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“in  this  study,  the  resolution  of  a  nonlinear  adjustment  model  is  addressed 
through  an  isomorphic  geometrical  setup  with  tensor  structure  and  notation, 
represented  by  a  u-dimensional  "model  surface""* embedded  in  a  flat  n-dimensional 
"observational  space’*.'  The  n  observations  correspond  to  the  observational-space 
coordinates  of  the  point  Q,  the  u  Initial  parameters  correspond  to  the  model- 
surface  coordinates  of  the  "initial”  point  P,  and  the  u  adjusted  parameters 
correspond  to  the  model-surface  coordinates  of  the  "least-squares"  point  P. 

The  least-squares  criterion  results  in  a  minimum-distance  property  implying 
that  the  vector  PQ  must  be  orthogonal  to  the  model  surface. 

The  geometrical  setup  leads  to  the  solution  of  "modified  normal  equations" 
characterized  by  a  positive-definite  matrix.  The  latter  contains  the  second 
order  and,  optionally,  the  third-order  partial  derivatives  of  the  observables 
with  respect  to  the  parameters.  This  approach  significantly  shortens  the 
convergence  process  as  compared  to  the  standard  (linearized)  approach.  The 
variance-covariance  and  the  weight  matrices  of  adjusted  quantities  correspond 
respectively  to  the  associated  metric  tensors  and  to  the  metric  tensors  at  P, 
formulated  in  the  pertinent  Rlemannlan  manifolds. 
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1.  INTRODDCTION 


The  parametric  adjustment  model  expresses  each  of  the  observables  in  terms 
of  parameters,  where  the  structure  Jinking  the  two  groups  of  varlabies  is,  in 
general,  nonlinear.  The  number  of  observables  is  denoted  by  n  and  the  number  of 
parameters  by  u.  where  n  must  be  greater  than  u  for  an  adjustment  to  take  place. 
The  adjustment  model  is  written  as 

=  F(X^)  . 

where  and  are  the  sets  (column  vectors)  of  adjusted  observations  and 
adjusted  parameters,  respectively.  This  study  addresses  the  resolution  of  a 
nonlinear  model  through  an  isomorphic  geometrical  setup  with  tensor  structure 
and  notation.  Such  efforts  date  back  to  [Blaha,  1984],  which  treats  a  linear  or 
linearized  adjustment  model.  An  initial  analysis  of  a  nonlinear  model  can  be 
found  in  [Blaha,  1987]. 

In  a  standard  adjustment  approach,  a  nonlinear  adjustment  model  Is  subject 
to  the  Taylor  series  expansion  based  on  an  initial  set  of  parametric  values,  X° . 
The  terms  in  the  second  and  higher  powers  of  the  parametric  corrections  are 
neglected,  resulting  in  the  familiar  (linearized)  observation  equations.  In 
matrix  notation,  the  latter  are  expressed  by 

V  =  AX  ^  L  , 

where  A  is  the  design  matrix,  X’X^-X°  is  the  column  vector  of  parametric 
corrections.  V=L^-L^  is  the  column  vector  of  residuals,  and  L=L°-L^  is  the 
column  vector  of  constant  terms,  with  L°  =  F(x'’)  containing  the  values  of 
observables  consistent  with  the  initial  set  of  parameters  and  containing  the 
actual  observations.  The  linearized  model  is  subjected  to  the  least-squares 
criterion 

P  V  =  minimum  , 

where  P  is  the  weight  matrix  of  observations.  This  criterion  leads  to  the 
formation  of  the  familiar  normal  equations. 

If  the  original  adjustment  model  Is  nonlinear,  the  resolution  of  the 
linearized  model  does  not  yield  the  final  answers.  The  process  is  usually 
repeated  with  new,  updated  parameters  and  the  corresponding  changes  in  A  and  L. 
However,  the  variance-covariance  matrix  of  observations.  Z,  as  well  as  the 


weight  matrix  P,  adopted  as  P=£  are  Independent  of  the  parameters.  Thus,  the 

T 

matrix  of  normal  equations.  N  =  A  PA,  changes  only  due  to  A,  and  the  column 

T 

vector  representing  the  right-hand  side  of  normal  equations,  U“-A  PL,  changes 
only  due  to  A  and  L.  The  computation  of  the  updated  parametric  values  through  a 
new  X  requires  the  formation  and  the  Inversion  of  a  new  N  In  each  Iteration,  or 
a  mathematically  equivalent  procedure.  When  X  becomes  sufficiently  close  to 
zero  the  iterative  process  is  terminated.  As  its  by-product,  the  latest  matrix 
N  ^  is  adopted  as  the  variance-covariance  matrix  of  adjusted  parameters. 

The  functional  relationship  between  the  observables  and  the  parameters 
lends  itself  to  a  geometrical  interpretation  and  treatment  involving  spaces  and 
surfaces  generalized  to  higher  dimensions.  In  particular,  the  parametric 
adjustment  model  symbolized  by 

X  =•  X  (u  )  .  r  =  1,2 . n  ,  a  =  1,2, ...,u 

can  be  linked  to  the  Gauss  form  of  a  surface  in  relation  to  the  surrounding 
r  Cl 

space,  where  x  are  the  apace  coordinates  and  u  are  the  surface  coordinates. 

The  Gauss  form  of  a  two-dimensional  surface  (u»2)  embedded  in  a  three- 
dimensional  flat  space  (n=3)  is  described,  together  with  two  other  forms,  in 
Chapter  6  of  [Hotlne,  1969].  In  [Blaha,  1984],  both  the  n-dimensional 
"observational"  space  and  the  u-dlmensional  "model"  surface  were  considered 
flat.  The  latter  was  thus  in  reality  a  hyperplane.  Although  the  model  surface 
is  now  Intrinsically  a  curved  space,  the  surrounding  space  will  be  seen  to  be 
again  flat. 

In  general,  a  flat  space  can  be  described  via  Cartesian  coordinates.  In  a 
Cartesian  coordinate  system,  a  given  point  is  depicted  by  a  set  of  coordinates 

r 

X  ,  r=l,2,...,n,  which  can  be  interpreted  as  its  position  vector  expressed  by 

r  r 

contravariant  components,  p  =x  .  Considered  as  a  tensor,  and  thus  as  a  point 

P 

function,  p  is  associated  with  this  given  point  (and  not,  for  example,  with  the 

origin).  Accordingly,  if  the  space  coordinates  are  Cartesian,  the  above  model 
r  r  Of 

equation  x  =x  (u  )  will  express  a  family  of  position  vectors  associated  with  the 

model  surface.  Any  of  these  position  vectors  could  be  freely  parallel- 

transported  to  any  location  in  space,  e.g.  to  the  surface  point  P  described  by 
r  r  d 

x^-x  (u^),  and  could  eventually  give  rise  to  tensor  equations  there. 

Although  derived  in  Cartesian  coordinates,  such  equations  would  be  valid  In  any 
coordinates  applicable  to  a  flat  space. 


2 


In  an  Important  extension  of  the  above  discussion,  we  state  that  if 

g  =  constant  . 
sr 

r  ot 

where  g^^  is  the  space  metric  tensor,  x  (u  )  can  again  be  interpreted  as  a 
family  of  position  vectors  associated  with  the  model  surface  embedded  in  a  flat 
space.  Thus,  although  the  space  coordinate  system  cannot  be  general  for  such  an 
interpretation  to  hold  true,  it  is  not  required  to  be  Cartesian.  The  position 
vectors  can  again  be  parallel-transported  to  the  point  P  or  any  other  location 
without  changes  in  their  components  (contravariant  as  well  as  covariant),  which 
has  been  referred  to  as  "free"  parallel  transport.  This  stems  from  the  fact 
that  such  changes  are  expressible  in  terms  of  the  Chrlstoffel  symbols  formed 
through  partial  derivatives  of  the  space  metric  tensor. 

Among  ail  possible  space  coordinate  systems,  only  those  characterized  by  a 

constant  metric  tensor  will  be  relevant  to  our  development.  Accordingly,  x^iu**) 

will  be  Interpreted  as  a  family  of  position  vectors  associated  with  the  model 

surface,  any  of  which  can  be  freely  parallel-transported  to  a  chosen  location. 

The  foregoing  is  meant  to  elucidate  the  isomorphism  of  adjustments  and  geometry. 

At  the  same  time,  it  has  led  us  to  appreciate  the  qualitative  difference  between 

r  r  Ot 

the  two  contexts.  In  the  first  context,  the  model  equation  x  =x  (u  )  provides 
restrictions  on  adjusted  observations  of  a  general  kind,  where  ail  the 
quantities  and  relationships  can  be  completely  void  of  any  geometrical  meaning. 
In  the  second  context,  it  provides  restrictions  on  position  vectors  in  space, 
which  is  a  purely  geometrical  phenomenon.  However,  studying  one  scientific 
discipline  In  terms  of  another,  even  seemingly  disparate,  may  lead  to  unforeseen 
benefits  and  a  deeper  urderstandlng  of  both. 
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2.  GEONETKICAL  SBTOP 


J* 

In  denoting  the  n  observables  by  x  ,  r»l,2 . n,  and  the  u  unknown 

parameters  by  u**,  a=l,2,...,u,  wt  can  represent  a  nonlinear  parametric 
adjustment  model  by 

f  a.  r  .r  .  a  *  *  *  P 

X  =x(u)=x  +A  Au  +(l/2)naAu  Au 

o  <x  aip 

y  (1/6)  Au®  Au^  Au’'  +  ...  ,  (la) 

apr 

,  ct  a  <x  , , .  V 

Au  =■  u  -  u  ,  (lb) 

o 

where  u®  represents  an  Initial  set  of  parameters  and  x^=x^(u®)  represents 

O  0  0 

the  observables  consistent  with  this  set.  Throughout  this  study,  the  lower-case 
Roman  Indices  range  from  1  to  n,  and  the  lower-case  Greek  Indices  range  from  1 
to  u.  Tensor  symbolism  Implies  the  summation  convention  over  the  dummy 
(repeating)  Indices. 

In  the  geometrical  context,  the  first  equality  in  (la)  represents  the  Gauss 

form  of  a  u-dlmenslonal  surface  embedded  In  an  n-dimenslonal  space.  The  surface 

Is  endowed  with  the  coordinate  system  {u®},  a=l,2 . u  and  is  referred  to  as 

the  model  surface,  and  the  space  Is  endowed  with  the  coordinate  system  (x 

r=l,2 . n.  and  Is  referred  to  as  the  observational  space.  The  second  equality 

j.'’  (la)  is  the  Taylor  series  expansion  of  x*'  from  the  "Initial"  point  P  lying  in 

the  model  surface,  whose  model-surface  coordinates  are  u®  and  whose 

0 

observational -space  coordinates  are  x^.  The  notation  Identifying  the  partial 

derivatives  at  P,  such  as  8x^/0u®aA*^,  3  ^x''/3u®au^aO^. ,  etc.,  Is  self- 

a  otp 

evident.  Th"  actual  observations  can  be  thought  of  as  describing  the  point  Q  in 
the  observational  space,  which,  due  to  measuring  errors,  cfo"r.  not  lie  In  the 
known  model  surface.  The  task  at  hand  consists  in  determining,  from  the 
observed  point  Q,  a  model -surface  point  satisfying  the  least-squares  criterion. 

In  the  adjustment  context,  the  variance-covariance  and  the  weight  matrices 
of  observations  depend  on  the  quality  of  measurements.  They  are  independent  of 
the  adjustment  model,  of  the  initial  set  of  parameters,  of  the  outcome  of 
obser''atlon3 ,  etc.  Thus,  for  a  given  observational  mode  they  are  constant.  In 
the  "traditional"  identification  of  (Blaha,  1984],  variance-covariance  .matrices 
correspond  to  associated  metric  tensors,  and  weight  matrices  correspond  to 
metric  tensors.  Accordingly,  we  represent  the  variance-covariance  matrix  of 
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observations  by  the  observational-space  associated  metric  tensor  g  ,  and  the 
weight  matrix  of  observations  by  the  observational-space  metric  tensor  ■  and 
state  that  both  tensors  are  independent  of  the  form  of  the  model  surface,  of  the 
initial  point  P,  of  the  observed  point  Q,  etc.,  leading  to  the  simplification 

g  =  constant  .  (2) 


One  could  also  attribute  the  tensors  g  and  g^^  to  the  point  Q  and  state  that 
the  geometrical  setup  must  account  for  Q  located  anywhere  in  the  observational 
space.  In  turn,  (2)  implies  that  the  observational  space  must  be  flat. 

-  r* 

If  the  set  (5x  denotes  the  coordinate  differences  between  the  observed 
point  Q  and  the  desired  model -surface  point  denoted  P.  it  corresponds  to  the 
negative  residuals,  and  the  least-squares  criterion  corresponds  to 

<53^^  "  <5x^  g  <5x^  =  minimum  .  (3) 

®  g  r 

Since  (2)  allows  us  to  make  appeal  to  Cartesian  coordinates,  we  can  consider  the 

"  r  r 

quantities  such  as  (5x  ,  Ax  ,  etc.,  whether  Infinitesimal  or  finite,  to  be  sets 
of  coordinate  differences  as  well  as  contravar iant  vectors  in  the  observational 
space,  where  they  can  be  freely  parallel-transported  to  any  location.  .Moreover, 
the  quadratic  form  (3)  represents  the  square  of  the  distance  between  Q  and  P. 
Therefore,  the  desired  " least -squares"  point  P  must  be  the  foot-point  of  the 
straight  line  dropped  or''hogonal  ly  from  Q  onto  the  model  surface.  This  property 
is  schematically  Illustrated  in  Fig.  1.  where  the  vector  PQ.  l.e.,  the 
geometr  ic.il  object  symbolized  by  <Jx.  is  drawn  perpendicular  to  the  model 
surface  at  P.  We  note  that  if  any  other  adjustment  condition  were  used  in 
lieu  of  the  least  -squares  criterion,  the  minimum-distance  property  (3)  would  not 
exist  and  the  geometric-tensorlal  treatment  of  the  adjustment  th»r.fy  would 
probably  be  much  more  complex  if  not  impossible. 

The  model -surface  coordinates  of  the  least-squares  point  P,  denoted  u^, 

correspond  to  the  adjusted  parameters,  and  Its  observational -space  coordinates, 

denoted  x ,  correspond  to  the  adjusted  observations.  (Although  this  designation 

(X  r 

Is  more  restrictive  than  the  notation  u  and  x  used  In  la.b.  It  need  not  cause 
any  confusion.)  Furthermore,  the  surface  associated  metric  tensor  at  P 
corresponds  to  the  variance-covariance  matrix  of  adjusted  parameters,  and  the 
necessary"  version  of  the  space  associated  metric  tensor  at  P  corresponds  to 
the  var lar.cc-covar lance  matrix  of  adjusted  observations.  Such 
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observational 


space 
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nodej 


surface 


(Ax'  )  .  (Au) 


Ax '  .  Au 


model  plane 


Pig.  1 

Symbolic  representation  of  a  nonlinear  parametric  least-squares  adjustment. 

A  one-dlmens lonal  model  surfare  Is  embedded  In  a  two-dimensional  observational 
space;  the  parentheses  pertain  to  the  solution  of  the  model  linearized  at  P. 
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correspondences,  as  well  as  the  definition  of  "necessary"  tensors,  were 
presented  in  Sections  2.2  and  4.2  of  [Blaha,  1984].  In  the  same  vein,  the 
surface  metric  tensor  at  P  corresponds  to  the  weight  matrix  of  adjusted 
parameters,  and  the  "necessary"  version  of  the  space  metric  tensor  at  P 
corresponds  to  the  singular  weight  matrix  of  adjusted  observations. 

These  notions  confirm  the  variance-covariance  propagation  law  and  Introduce 
a  weight  propagation  law",  and  they  can  be  extended  to  functions  of  adjusted 
parameters  as  well  as  to  functions  of  adjusted  observations.  In  [Blaha,  1984], 
however,  they  were  demonstrated  only  for  linear  models.  In  nonlinear  models  the 
above  tensors,  expressible  by  means  of  orthonormal  vectors  tangent  to  the  model 
surface  at  P.  could  be  said  to  describe  "tangential",  or  first-order, 
properties  of  least-squares  estimates. 
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3.  SOLDTION 


3 . 1  General  Development 

The  adjustment  theory  assumes  that  the  known  values  uf  the  Initial  set  of 
parameters,  represented  here  by  u*.  are  close  to  the  final  values  represented 
by  u^.  Furthermore,  the  values  of  the  residuals  are  also  assumed  to  be  small. 
For  the  needs  of  the  current  development,  we  stipulate  that  quantities  will  be 
considered  "small",  or  first-order.  If  they  behave  like  differential  quantities, 
in  the  sense  that  a  product  of  two  first-order  quantities  will  result  In  a 
quantity  smaller  by  an  order  of  magnitude.  The  latter  will  be  referred  to  as  a 
second-order  quantity.  In  the  geometrical  context  the  assumption  about  u*  is 
maintained,  implying  that  the  model-surface  coordinate  differences  between  P 
and  P,  in  (lb)  symbolized  by  Au**,  are  small.  The  same  assumption  applies  also 
for  the  observational-space  coordinate  differences  between  these  two  points, 

A  r"  f  I  A\ 

Ax  =  X  -  X  .  ( 4 ) 

o 

However,  the  values  In  <5x''  =  x^-x^,  corresponding  to  the  negative  residuals, 
are  subject  to  no  such  restriction.  This  is  the  first  distinction  between  the 
standard  and  the  geometrical  approaches  to  nonlinear  models. 

The  second  distinction,  which  Is  the  cornerstone  of  our  geometrical 
approach.  Is  a  consequence  of  (3).  Before  drawing  conclusions  from  this 
minimum-distance  criterion,  we  review  the  notion  of  the  design  tensor, 
corresponding  in  adjustment  calculus  to  the  design  matrix.  In  (Blaha,  1984]  and 
in  Section  2.2  of  [Blaha,  1987],  this  tensor  was  identified  with  the  set  of 
partial  derivatives  9x^/5u*  at  the  initial  point  P,  and  was  presented  in  a 

vectorized  form.  Thus.  In  (la)  Is  the  design  tensor  at  P.  which  can  be 

<x 

expressed  by  means  of  a  set  of  orthonormal  vectors  spanning  the  hyperplane 
tangent  to  the  model  surface  at  P: 

a’’  3  ax^'/au®  =  e''  4  +  j*'  j  +  ...  .  (5) 

a  a  a 

The  above  .hyperplane  is  referred  to  as  the  "model  plane",  and  In  Fig.l  it 

r  r 

is  symbolically  represented  by  a  straight  line.  The  notation  8  .  j  ,  ... 
Identifies  the  orthonormal  vectors  at  P  spanning  the  model  plane  by  their 
observational -space  contravar iant  components,  while  the  notation  8^.  j^.  ... 
identifies  them  by  their  model -surface  covariant  components.  In  general,  when 
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Identifying  individual  vectors  as  geonetrical  objects,  we  remove  the  indices. 

The  double  Identification,  such  as  Ax'  and  Au  written  side-by-slde  in  Pig.  1, 
depicts  one  and  the  same  vector  exiu-essed  in  two  coordinate  systems.  The 
notions  of  model  surface  and  model  plane  are  further  discussed  In  Appendix  B, 
which  also  describes  in  detail  the  role  of  the  design  tensor  in  either  manifold. 


The  design  tensor  (5)  serves  in  formulating  the  model -surface  metric  tensor 
at  the  initial  point  P  as 


3a 


•^3  '^a 


=  A.  g  A 
3  ^sr  a 


(6) 


.‘delations  similar  to  (5)  and  (6)  appear  in  Chapter  6  of  [Hotlne,  1969J,  where 
they  are  derived  for  n=3  and  u=2.  Equation  (5)  provided  an  early  impetus  In 
treating  the  adjustment  theory  via  geometry.  I^'s  motivating  power  stems  from 
the  fact  that  its  left-hand  side  describes  the  algebraic  structure  of  the  design 
matri.x  in  tne  parametric  adjustment  model,  and  its  right-hand  side  gives  this 
structure  a  clearcut  geometrical  meaning. 


r 

In  analogy  to  A^,  the  design  tensor  at  P  can  be  formulated  as 

.p  «p  .-  -p- 

A  =  e  fi  *  j  J  *  . 

a  a  a 

where  the  overbars  identify  a  set  of  orthonormal  vectors  spanning  the  hyperplane 
tangent  to  the  model  surface  at  P.  Due  to  the  orthogonality  of  (Sx  and  the 
model  surface  as  dictated  by  (3),  at  P  we  must  have 


^3  ®sr  ° 

S  —  —  p 

Since  A  is  a  mixed  space-surface  tensor  at  P.  and  <5x  Is  a  contravariart 

3  _ 

space  vector  at  P,  the  above  Is  a  tensor  equation  at  P.  This  equation 
appeared  in  a  similar  form  *n  Section  2.3  of  [Blaha,  1987J.  The  orthogonality 
condition  (7)  embodies  the  above-mentioned  second  distinction,  in  that  it  is  a 
precise  relation  whose  active  exploitation  does  not  have  an  equivalent  in  the 
standard  adjustment  approach.  As  has  been  indicated  in  the  Introduction,  the 
latter  proceeds  by  iterations,  where  a  set  of  values  approximating  Au**  is 
refined  at  each  step.  When  the  changes  are  sufficiently  small,  the  iterative 
process  terminates  and  (7)  is  satisfied  as  a  by-product.  However,  upon  seeking 
to  fulfill  this  condition  at  each  individual  step,  the  standard  approach  can  be 
modified  and  the  iterative  process  can  be  significantly  shortened. 
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since,  in  the  geometrical  context,  P  la  assumed  to  be  close  to  P,  a  good 
approximation  to  Au**  is  offered  by  (Au**),  where  the  parenthesis  notation 
represents  the  solution  of  the  linearized  model  described  by  (la.b)  with 

0^-  =  0,  =0 .  This  approximation  is,  in  fact,  the  initial  iteration 

of  the  standard  approach.  We  note  that  if  the  values  in  Ax  are  so  large  as 
to  approach  finite  quantities,  a  better  approximation  to  Au”'  can  be  obtained  in 
a  different,  more  general  form  presented  in  the  next  section.  In  the  current 
discussion  the  linearized  model  corresponds  to  the  situation  where  the  model 
plane  replaces  the  model  surface.  We  retain  {u**}  as  the  model-plane  coordinate 
system,  and  observe  that  (6)  is  the  plane's  constant  metric  tensor.  In  analogy 
to  a  statement  made  earlier  in  conjunction  with  the  constant  space  metric 
tensor,  we  can  consider  the  quantities  such  as  (Au®),  Au® ,  etc.,  whether 
Infinitesimal  or  finite,  to  be  sets  of  model-plane  coordinate  differences  as 
well  as  contravar lant  vectors  in  the  model  plane,  where,  if  need  be.  they  could 
be  freely  parallel-transported  to  any  location. 

p 

The  model-plane  covariant  components  of  the  orthogonal  projection  of  <5x 
onto  the  model  plane  are  obtained  with  the  aid  of  the  design  tensor  (5); 

hr 

where  the  contravar iant  vector  Ax^.  known  a  priori,  can  be  decomposed  as 


r  r  r  r  r 

Ax  s  Xq  -  =  (Ax'  )  +  (Ax"  )  . 

In  terms  of  Fig.  1,  the  double  identification  (Ax')  and  (Au)  designates  one  and 
the  same  vector  PPj •  Based  on  (8),  we  have 

(Au®)  =  a®^  (^u„)  ,  (9) 

P 


where  a  ^  is  the  mode  1 -surface  associated  metric  tensor 
model -plane  associated  metric  tensor,  obtained  from  a^^ 


a®^ 


at  P,  or  the  constant 
according  to 


If  (Au®)  la  considered  to  be  a  set  of  model-surface  rather  than  model -plane 
coordinate  differences,  it  gives  rise  to  a  surface  point  which,  in  Fig.  1,  is 
denoted  (P).  In  consulting  the  figure,  we  note  that  in  the  observational  space 
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this  point  Is  arrived  at  by  neans  of  the  vector  (Ax),  where 


(Ax^)  =  (Ax'*')  +  (v*')  . 

Topics  similar  to  the  above  are  discussed  In  Appendix  C. 


If  the  adjustment  model  were  linear,  the  model  surface  would  reduce  to  the 
model  plane  and  (Au*)  would  represent  the  final  parametric  corrections.  The 
design  tensor  (5)  would  serve  to  express 

(Ax'^')  =  (Au“)  . 

r  r 

which  would  yield  the  adjusted  observations  as  x^+(Ax'  ).  The  negative 
residuals  would  be  found  from 


(Ax"'')  =  -  (Ax''')  . 


The  variance-covariance  matrix  of  adjusted  parameters  would  correspond  to  a**^, 
while  the  var lance-covariance  matrix  of  adjusted  observations  would  correspond 
to  the  necessary  associated  metric  tensor: 


g-"  - 


>  J"  J* 


(10a) 


The  weight  matrix  of  adjusted  parameters  would  correspond  to  a^^  and  the  weight 
matrix  of  adjusted  observations  would  correspond  to  the  necessary  metric  tensor: 


■  ®sl  Sjr  • 

Furthermore,  the  variance-covariance  and  the  weight  matrices  of  residuals  would 
correspond  respectively  to 


„r3  rs  .rs 
S  -  g  -  g 


g"  =  g  ~  g' 
sr  ®sr  sr 


The  above  relations  were  presented  In  [Blaha,  1984]. 


(lla.b) 


In  the  current  nonlinear  model,  (Au®)  represents  the  first  approximation  to 
Au*.  We  remark  that  In  a  two-  or  higher-dimensional  model  plane  the  vectors 
(Au)  and  Au  are  not  colllnear  in  general.  The  set  Au*.  when  considered  In  the 
role  of  model-surface  coordinate  differences,  gives  rise  to  the  desired  least- 
squares  point  P  (see  Fig.  1).  In  the  observational  space  this  point  Is 
obtained  by  means  of  the  vector  Ax,  where 


r  r  r 

Ax‘  “  Ax'  +  V  . 


(12a) 


1  1 


ProM  (la.b)  and  (4),  we  express  the  right-hand  side  of  (12a)  by 

Ax*''  »  Au“  ,  (12b) 

a 

v*"  =•  (1/2)  Au*  Au^  +  (1/6)  Au*  Au^  Au^  +  ...  .  (12c) 

ap  apr 

This  outcome  is  also  discussed  In  Appendix  C. 

Whether  the  least-squares  point  P  should  be  described  by  the  model- 
surface  coordinates  or  by  the  observational-space  coordinates  via  (12a-c),  the 
values  Au^  are  unknown  and  unobtainable  by  closed  formulas  sucl.  as  (8)  and  (9) 
applicable  to  the  linear  model.  In  view  of  formulating  a  useful  relationship 
for  Au*.  we  first  observe  that 

A*”  «  a'*  H-  AA*'  ,  (13a) 

a  a  a 

Aa"^  =  o'".  Au^  +  (1/2)  o'*.  Au^  Au’'  +  ...  .  (13b) 

a  ap  apj 

r  S  V 

where  (13b)  stems  from  the  fact  that  at  P  we  have  3A  /9u  =0 

.  a  ap  _ 

3^A''/3u^au^  =  o'^.  ,  etc.  Equations  (13a, b)  express  the  design  tensor  at  P 
from  Its  counterpart  at  P.  Further,  Fig.  1  shows  that  Ax  can  be  expressed  as 

Ax*^  =  Ax*^  -  Ax*^  ,  (14) 

where  Ax'*  Is  known  but  Ax'^,  determined  via  Au®  in  (12a-c),  is  unknown.  This  is 

-  r  r 

a  tensor  equation  at  the  least-squares  point  P,  to  which  Ax  and  Ax  can  be 
freely  parallel-transported  from  the  initial  point  P.  The  tensors  A^  and 
Ax'*  at  P  are  needed  in  view  of  (7). 

The  approximation  set  (Au*)  will  be  corrected  by  A(Au*)  to  yield  the 
desired  set  Au* : 

Au*  =  (Au**)  +  A(Au®)  .  (15) 

The  unknown  set  A(Au“)  will  be  determined  through  the  precise  relation  (7).  We 
stipulate  that  the  level  of  accuracy  In  (7)  should  be  of  the  third  order,  in  the 
sense  that  the  neglected  terms  can  only  be  of  the  fourth  and  higher  orders.  By 
a  hypothesis  common  to  both  the  standard  and  the  geometrical  approaches,  the  set 
Au®  contains  small  quantities  and  the  set  (Au®)  represents  the  main  contribution 
of  Au®.  Thus,  (Au®)  and  Au®  are  assumed  to  contain  first-order  quantities, 
while  the  set  A(Au®)  is  assumed  to  contain  second-order  quantities.  Finite 
quantities  contracted  with  (Au®)(Au^)  are  considered  to  be  of  the  second  order. 


similar  to  finite  quantities  contracted  with  AiAu**).  Fourth-order  quantities  to 
be  neglected  are  of  three  kinds;  o[(Au°')  (Au^)  (Au^)  (Au*^)  ] ,  o[(Au*)(Au^)  A(Au^)), 
and  ofAlAu**)  A(Au^)].  Since  the  last  kind  Involves  quadratic  terms  in  A(Au*) , 
the  reason  for  carrying  out  the  solution  to  the  third  order  becomes  clearer. 


In  considering  the  known  quantities  (Au  ),  similar  to  (13a, b)  we  have 

(16) 


(A^)  =  ^  0^^  (Au^)  +  (1/2)  (Au^)  (Au’')  +  ... 


-r  r 

We  now  form  A  based  on  the  outcome  for  (A  ) : 

a  01 

A^  -  (A^)  ^  A(A^)  . 
a  a  a 


(17a) 


r  ”  r  CL 

The  set  A(A  )  can  be  obtained  as  A  in  (13a, b)  with  Au  substituted  for  by 
a  «  ' 

the  right-hand  side  of  (15),  minus  (A**)  from  (16).  This  yields 


A(A^)  =  0^^  A(Au^)  >  (Au’')  A(Au^)  +  ...  . 

P 

where  use  has  been  made  of  the  symmetry  of  .  In  a,  i3.  and  r  through 

ap7 

A{Au^)  (Au’')  *■  (Au^)  A(Au’')  =  2  (Au’')  A(Au^)  . 


(17b) 


a/Sr  '  '  '  '  a0r 

The  dots  in  {17b)  represent  o[A(Au’^)  A(Au^)  ] +o[  (Au*^ )  (Au’' )  A(Au^)l  +  .,.,  which 
are  respectively  fourth-,  fourth-,  and  higher-order  quantities. 


Furthermore,  similar  to  (12a-c)  we  write 

(Ax^)  =  (S  (Au®)  +  (1/2)  o''.  (Au“)  (Au^) 
a  ixp 

*  (1/6)  (Au®)  (Au^)  (Au’')  +  ...  .  (18) 

apT 

r  r 

In  analogy  to  the  above,  we  form  Ax  based  on  the  outcome  for  (Ax  ): 

Ax'"  =  (Ax'")  +  A(Ax'')  .  (19a) 

r  r  Cl 

The  set  A(Ax  )  can  be  obtained  as  Ax  in  (12a-c)  with  Au  given  by  (15),  minus 

(Ax'^)  given  by  (18).  This  yields 

A(Ax'’)  =  a''  A(Au®)  +  o'".  (Au^)  A(Au“)  <-...,  (19b) 

a  ap 

P 

where  we  have  employed  the  symmetry  of  0  „  in  a  and  0  through  the  identity 

ap 

0^^  A(Au®)  (Au^)  +  0^^  (Au®)  A(Au^)  -  2  0^^  (Au^)  A(Au®)  . 
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However,  due  to  (16),  equation  (19b)  simplifies  to 

A(Ax'‘)  =  (A^)  A(Au“)  (19c) 

a 

where  the  dots  represent  o[A(Au^)  ACAu**)  J +  o[(Au^)  (Au^)  (Au^)  A(Au®)J  +  ...,  which 
are  respectively  fourth-,  fifth-,  and  higher-order  quantities. 

—r 

Finally,  due  to  (14)  and  (19a),  we  can  express  Ax  as 

Ax'*  -  (Ax'*)  -  A(Ax‘*)  ,  (20a) 

where 

(Ax'*)  =  Ax'*  -  (Ax'*)  .  (20b) 

In  terms  of  {17a)  and  (20a).  the  precise  relation  (7)  reads 
[(A^)  ^  A{A®)]  [(Ax'*)  -  A(Ax'*)]  =  0  . 

or 

hr  ®3r 

-  (A®)  gg^  A(Ax'')  -  A(A®)  A(Ax’‘)  =  0  .  (21) 

The  second,  third,  and  fourth  terms  on  the  left-hand  side  of  (21)  are  denoted 
respectively  as  a,  b,  and  c,  and  will  now  be  evaluated. 

In  considering  (17b)  we  obtain 

a  =  A(Au^)  (Ax'*)  «®^^  (Au*^)  A(Au^)  g^^  (Ax'*)  +  ...  . 

In  analogy  to  the  statement  that  followed  (17b),  the  dots  represent  fourth-  and 
higher-order  quantities.  However,  (Ax'*)  is  excluded  from  influencing  this 
order.  Further,  in  consulting  (19c)  we  write 

b  =  (Ap  gg^  (A^)  AlAu'*)  >  ...  . 

In  analogy  to  the  statement  that  followed  (19c),  the  dots  again  represent 
fourth-  and  higher-order  quantities.  Finally,  in  consulting  both  (17b)  and 
(19c)  we  conclude  that 

c  =  o[A(Au^)  AlAu**)]  +  ...  , 

which  Is  composed  entirely  of  fourth-  and  higher-order  quantities. 


Upon  collecting  the  results,  (21)  yields 

(a|)  g,^  (dx'')  =■  (A®)  gg^  (A^)  A(Au*)  -  A(Au“)  g^^  (fix'’) 

-  (Au'')  A(Au®)  gg^  (Ax'")  +  ...  .  (22) 

where  the  dots  represent  fourth-  and  hlgher-ordcr  quantities  regardless  of 
(dx**).  Equation  (22)  Is  equivalent  to  the  precise  relation  (7).  If  the  terms 
symbolized  by  dots  are  omitted,  the  equivalence  is  valid  to  the  third  order.  In 
this  case,  the  desired  set  A(Au*)  is  obtained  by  solving  the  linear  system 

"V  -  *lar  ■  **6'  '■’’■s*  ■ 

where 

'=Sa>  ■  “e’  '<>  ■ 

All  the  quantities  In  parentheses  are  formed  with  (Au*).  The  relation 

{23a)  could  be  called  "modified  normal  equations",  and  the  matrix  form  of  the 

quantity  inside  the  brackets  could  be  called  "matrix  of  modified  normal 

equations".  The  second-  and  the  third-order  partial  derivatives  grouped 

respectively  In  of  and  are  evaluated  at  P  and  can  thus  be  stored  and 

pa  par 

treated  as  constant  should  (23a)  be  re-applied.  Their  "once  and  for  all" 
evaluation  is  useful  since  In  complex  models  their  re-computatlon  at  an  updated 
point  would  represent  a  serious  drawback.  We  remark  that  if  (<5x  )  is  small,  in 
the  sense  it  contains  first-order  quantities,  the  third  term  Inside  the  brackets 
of  (23a)  can  be  left  out  since  it  would  give  rise  to  fourth-order  quantities. 


3 . 2  Initial  Iteration  with  Large  Residuals 


The  development  in  the  previous  section  was  based  on  the  set  (Au*)  as 
obtained  in  the  initial  iteration  of  the  standard  adjustment  approach.  In 
[Blaha,  1987],  such  a  common  Iteration  was  called  "zero-th  iteration".  This 
concept  is  abandoned  since  the  two  approaches  can  now  be  considered  distinct 
already  at  the  initial  stage.  In  particular,  we  develop  a  general  case  of  the 
geometrical  approach,  where  already  the  first,  or  Initial,  iteration  Is  based  on 
the  exact  relation  (7)  valid  at  the  least-squares  point  P. 


15 


In  consulting  the  previous  section,  we  present  the  barred  tensors  In  (7): 

A®  =  A®  +  Au’'  +  (1/2)  Au^  Au^  +  . .  .  ,  (24a) 

<5x^  -  Ax**  -  A^  Au“  -  (1/2)  0*“  Au“  Au®  -  .  . .  .  (24b) 

at  ae  '  ' 

To  avoid  expressions  nonlinear  In  Au*.  we  Introduce  approximations  In  (7) 
whereby  all  terms  containing  the  second  and  higher  powers  In  Au®  will  be  left 
out,  l.e.,  o(Au®  Au^),  o(Au®  Au^  Au^),  etc.,  will  be  neglected.  This  means  that 
the  third  and  further  terms  on  the  right-hand  sides  of  (24a, b)  will  be 
neglected,  as  well  as  the  term  -Q®  Au^g  A*^Au*  participating  in  (7). 

fi  7^  sf  (X 

Equation  (7)  is  thus  written  as 


hr  ~  ®sr  "a  ^  ‘^r  Sx''  *  .  .  .  =  0  ,  (25) 

where  the  dots  represent  the  terms  which  would  contain  the  second  and  higher 
powers  of  Au®,  and  where  no  assumptions  have  been  made  yet  about  the  size  of  the 
known  contravarlant  vector  Ax"^.  We  add  that  the  values  constituting  A^, 
and  g  are  assumed  finite  and  of  mutually  comparable  magnitudes. 


Upon  the  use  of  the  symbols  a^^  and  Ax^,  where 


0a 


A®  e  A** 

0  ®sr  a  ' 


Ax 


Son 

sr 


equation  (25)  becomes 


[a  -  Ax  Q®  <■...]  Au®  =  A®  Ax 
j3a  3  |3a  /3  s 


(26a, b) 


(27) 


The  first  geometrical  Iteration  Is  characterized  by  the  neglect  of  the  second 
and  higher  powers  of  Au®  in  (27),  In  which  case  the  solution  of  Au®  Is  denoted 
temporarily  as  (Au®); 


"  ^^3  ==  ^0  ^*s  •  ^27') 

In  terms  of  Fig.  1,  this  solution  would  lead  to  a  model-plane  point  lying 
somewhere  between  and  P^,  and  to  a  model-surface  point  lying  somewhere 
between  (P)  and  P.  If  we  further  neglect  -Ax^n®^Au®  and  denote  the 
corresponding  solution  by  (Au®),  equation  (27)  becomes 


“s«  ■  ''s  ■*’’3 


(27- I 


which  symbolizes  the  normal  equations  in  the  standard  adjustment  approach. 
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In  the  standard  approach,  the  Initial  values  of  parameters,  u*.  are 
assumed  to  be  close  to  the  final  values,  u  .  In  our  current  notation  and 
terminology,  this  assumption,  which  we  shall  retain,  is  equivalent  to  stating 
that  Au*^  contains  first-order  quantities.  Since  only  the  second-  and  higher- 
order  quantities  are  ignored  in  (27'),  {Au*}  in  the  first  geometrical  Iteration 
will  recover  the  main  contribution  of  Au®.  The  errors  in  {Au**}  will  thus  have 
the  magnitude  of  second-order  quantities,  consistent  with  the  property 
previously  attributed  to  the  initial  standard  iteration. 

As  the  comparison  of  (27")  with  (27')  reveals,  this  property  of  the  initial 

s  ct 

standard  iteration  is  indeed  valid,  provided  the  term  -Ax^n^^Au  neglected  in 
the  standard  approach  is  at  the  level  of  second-order  quantities,  i.e.,  at  the 

p 

level  of  the  terms  already  neglected  by  (27').  In  general,  this  is  true  if  i5x 

contains  first-order  quantities,  which  is  another  assumption  inherent  in  the 

standard  approach.  Accordingly,  under  the  two  standard  assumptions,  requiring 
ot  r 

that  both  Au  and  <5x  contain  first-order  quantities,  the  solution  of  (27")  is 
valid  to  the  first  order,  similar  in  this  respect  to  the  solution  of  (27').  The 
errors  in  both  (Au*)  and  (Au*)  are  then  second-order  quantities. 

p 

If  the  second  assumption  above  is  not  satisfied  and  Ax  contains  finite 
quantities,  (27')  will  still  hold  true  to  the  first  order  whereas  (27")  will 

Of 

not.  Although  the  quantities  in  (Au  ,  will  be  of  the  first  order,  their  errors 
will  In  general  attain  the  first  order  as  well,  due  to  the  neglect  of  the  first- 

3  ^ 

order  quantities  -Ax  fl.  Au  .  This  case  is  mainly  of  theoretical  interest. 

3  pa 

since,  in  practice,  it  corresponds  to  excessive  observational  errors  which 
should  be  eliminated  beforehand.  However,  we  might  also  consider  the  Instances 

p 

where  Ax  belongs  to  a  category  "in-between"  and  contains  quantities  that  are 
finite  but  small,  without  being  small  enough  to  warrant  the  attribute  "first 
order".  Here  the  solution  of  (27')  is  better  than  the  solution  of  (27"),  but 
not  by  an  order  of  magnitude.  We  conclude  that  in  general  the  solution  of  (27') 
approximates  Au*  as  well  as,  or  better  than,  the  solution  of  (27").  Better 

p 

results  are  obtained  when  the  values  in  Ax  have  the  magnitude  greater  than  that 
characterizing  the  first  order. 

p 

The  quantities  forming  Ax  and  their  link  to  (27')  and  (27")  can  be  further 
elaborated  on  as  follows.  Since  Au*  is  of  the  first  order  by  assumption,  and 
since  the  same  applies  to  (Au*)  and  (Au*).  equations  (27')  and  (27")  reveal  that 


the  tensor  A® fix  ,  which  figures  on  their  right-hand  sides,  is  also  of  the 
p  s 

first  order.  In  fact.  (27")  allows  this  tensor  to  be  written  as 
»ar  ■ 

The  above  expression  identifies  the  projection  of  the  vector  dx  onto  the  model 

plane,  as  previously  depicted  in  Fig.  1.  Since  (Au.)  is  of  the  first  order,  it 

P 

follows  that  if  the  magnitude  of  the  vector  <5x  is  significantly  greater  than 
that  characterizing  the  first  order,  this  vector  must  be  fairly  nearly 
orthogonal  to  the  model  surface.  This,  in  turn,  imposes  restrictions  on  the 
observational  errors. 

In  summary,  if  dx  is  of  the  first  order,  which  corresponds  to  the  presence 
of  small  observational  errors  as  encountered  in  practice,  the  first  geometrical 
iteration  (27')  has  the  same  level  of  accuracy  as  the  first  standard  Iteration 
(27").  This  category  Includes  dx^=0,  which,  in  turn.  Includes  the  case  of 
error-free  observations.  If  the  characteristics  of  dx  fall  into  the  "in- 
between"  category,  which  corresponds  to  the  presence  of  large  observational 
errors,  (27')  has  a  higher  level  of  accuracy  than  (27"),  although  not  by  an 
order  of  magnitude.  Here  the  errors  are  restricted  in  the  sense  that  dx  must  be 
fairly  nearly  orthogonal  to  the  model  surface.  Finally,  if  dx**  contains  finite 
values,  which  corresponds  to  the  presence  of  excessively  large  observational 
errors  (blunders)  and  is  thus  mainly  of  theoretical  interest,  (27')  is  still 
valid  to  the  first  order,  whereas  (27")  may  give  completely  false  results.  In 
this  case,  the  errors  are  restricted  in  the  sense  that  dx  must  be  very  nearly 
orthogonal  to  the  model  surface.  Since  (27')  gives  results  that  are  as  good  as 
those  of  (27")  in  most  cases,  and  better  in  the  remaining  cases,  the  use  of  this 
equation  as  the  first  geometrical  iteration  is  fully  Justified. 


3 . 3  Minimum-Distance  Property 
_2 

For  the  sake  of  theory,  we  express  ds  (the  square  of  the  distance 

2 

between  P  and  Q)  by  means  of  the  Taylor  series  based  on  ds  (the  square  of  the 
distance  between  P  and  Q).  The  latter  is  a  known  constant,  namely 

ds^  »  dx^  g  Sx^  , 
sr 
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while  the  former  is  a  function  of  Au  ,  expressed  as 

.0  — Q  _  P 

<Ss  »  5x  g  (5x  , 

sr 

~r  r 

where  6x  is  given  e.g.  in  (24b).  In  using  the  symmetry  of  g^^  and  0^^  in 

the  lower  indices,  we  deduce  that 

Ss^  =  Ss^  +  q  Au®  +  (1/2)  H  .  Au*  Au^  +  ...  , 

<x  ap 

where 

q  =  aAs^/au®  “  -  2  A®  g  Ax*^  , 
a  a  sr 

H  .  3  a^As^/au®au^  =  2  (A^  g  A^  -  Ax*  g  0^.) . 

a$  0  ^sr  a  ^sr  a$ 

The  second-order  partial  derivatives  grouped  in  are  assumed  continuous,  at 
least  in  the  neighborhood  of  the  point  P, 

A  similar  expression  can  be  written  for  a  model -surface  point  in  the 

2 

neighborhood  of  P,  whose  square  of  the  distance  from  Q,  denoted  As  ,  can  be 

-2 

expressed  in  terms  of  As  ; 

AS^  =  Ss^  *  q  Au“  +  (1/2)  H  .  Au®  Ad'®  ...  . 

(X  ap 

Here  Au*  Is  the  difference,  In  model-surface  coordinates,  between  the  point  in 

2 

question  and  the  point  P.  The  necessary  and  sufficient  conditions  for  A§ 
to  have  a  local  minimum  at  P  are 

q  =  -  2  A*  g  Ax*'  =  0  ,  (28) 

a  a  sr 


[H  „]  =  positive-definite  matrix  : 
ap 

the  brackets  indicate  the  matrix  equivalent  of  ft  where 

ap 

ft  -  2  (.A*  g  A^  -  Ax*  g  qI  )  . 
ap  p  sr  a  sr  ap 

The  latter  equation  can  be  written  as 


=  a^^  -  Ax^  . 


a.  »  Aa  g  A^  ,  Ax  »  g  Ax'^  . 
0a  0  ^sr  a  s  sr 


(31a,b) 
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According  to  (29),  the  above  quantity,  symmetric  In  the  Indices  a  and  should 
represent  a  positive-definite  matrix.  If  we  now  compare  the  right-hand  side  of 
(33)  with  the  quantity  Inside  the  brackets  of  (23a)  In  Its  final  update,  l.e., 

V  . 

In  the  last  Iteration  where  (a.  ),  (<Jx  ),  and  (Au  )  become  essentially  a-  , 

pa  s  pa 

<5x^.  and  Au^,  we  observe  that  within  the  required  accuracy  the  two  forms  are 
Identical.  It  Is  thus  confirmed  that  the  algorithm  represented  by  (23a)  leads 
to  a  minimum  provided  the  updated  matrix  of  modified  normal  equations  is 
positive-definite,  and  vice-versa.  Since  the  existence  of  a  minimum  Is  assumed 
throughout,  and  since  P  is  assumed  to  be  sufficiently  close  to  P,  the  above 
outcome  implies  that  the  matrix  of  modified  normal  equations  will  be  positive- 
definite  in  all  Iterations,  with  or  without  the  contribution  of  the  third-order 
partial  derivatives.  This,  In  turn,  Implies  that  it  can  be  inverted  by  the 
Choleskl  algorithm,  as  is  illustrated  In  the  example  presented  In  Appendix  A. 
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4.  DISCUSSION 


4 . 1  Basic  Considerations 

The  leading  subject  of  this  discussion  Is  equation  (23a).  applied  after  the 

r 

initial  Iteration.  We  begin  by  relating  (Sx  )  to  Sx  ,  where,  in  the 

adjustment  context,  the  latter  represents  the  negative  residuals.  From  (20a)  we 

write 

(dx*^)  -  <5x^  =  A(Ax^)  . 
which,  upon  using  (19c),  becomes 

(dx'')  =  dx'"  >  (a'*)  A(Au“)  +  ...  .  (34) 

a 

r 

Here  the  dots  again  represent  fourth-  and  higher-order  quantities.  Thus,  (Sx  ) 
differs  from  dx^  by  second-  and  higher-order  quantities. 

We  distinguish  three  basic  cases  with  regard  to  the  magnitude  of  the 

-r 

quantities  comprising  the  set  dx  : 

“  r  n 

1)  If  dx  is  of  the  second  order  or  smaller,  (dx  )  Is  In  general  of  the 

second  order  and  the  two  sets  (of  very  small  values)  may  be  quite  different. 

This  case  Includes  dx’'»0,  where  the  point  Q  lies  on  the  model  surface  as  would 
happen,  for  example,  with  error-free  observations. 

2)  If  dx^  is  of  the  first  order,  so  is  (dx'').  and  the  two  sets  are  similar 
in  values. 

3)  If  dx^  Is  large  (finite),  so  is  (dx*').  and  the  two  sets  are  nearly  equal 
in  values. 

In  the  first  case,  (dx  )  is  of  the  second  order  and  (23a)  simplifies  to 

s 

(a^^)  A(Au®)  =  (A^)  (dx_)  .  (35) 

This,  however,  is  precisely  the  formula  for  standard  adjustment,  corresponding 

to  normal  equations  formed  after  the  initial  Iteration.  In  particular,  the 

left-hand  side  is  A(Au„),  and  the  equation  corresponds  to  the  system  (8) 

P 

considered  at  the  point  (P).  We  conclude  that  in  this  theoretical  case,  which 
includes  the  category  of  error-free  observations,  the  standard  solution  is  valid 
to  the  third  order  and  hence  conforms  to  the  strategy  pursued  in  this  study. 
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The  second  case,  which  satisfies  the  standard  assumption  that  both  the 
parametric  corrections  and  the  residuals  are  small  (first-order),  is  certainly 
the  most  significant  in  practice.  As  was  Indicated  at  the  close  of  Section  3.1, 
equation  (23a)  becomes 

~  A(Au®)  =  (a|)  (dx^)  .  (36) 

In  comparing  it  with  the  standard  equation  of  the  form  (35),  we  notice  that  the 
latter  entails  the  error  (dx^)  fl^^A(Au*) ,  which  is  of  the  third  order.  The 
standard  formula  In  such  practical  cases  is  thus  valid  only  to  the  second  order. 
The  newly  developed  equation  (36)  is  superior  to  the  standard  normal  equations 
by  one  order  of  magnitude,  which  leads  to  a  commensurate  Increase  in  the  rate  of 
convergence.  In  fact,  depending  on  the  smallness  of  Au**,  and  thereby  on  the 
proximity  of  the  initial  point  P  to  the  least-squares  point  P,  the  solution  of 
( 36 )  may  be  final . 

Next,  we  examine  the  Improvement  in  the  solution  if  equation  (36)  is  re¬ 
applied.  We  further  assume  that  the  original  solution  A(Au*)  Is  valid  to  the 
third  order,  the  same  as  the  accuracy  of  the  formula  (36)  itself.  The  new. 
updated  set  (Au®)  is  thus  valid  to  the  third  order  (Instead  of  the  first  order 
as  assumed  initially),  and  the  new  A(Au*)  contains  fourth-order  quantities.  At 
this  stage  (36)  reveals  that  fifth-order  quantities  are  retained,  while  slxth- 
and  higher-order  quantities  are  ignored.  The  new  solution  will  thus  be  valid  to 
the  fifth  order.  By  contrast,  the  second  iteration  in  the  standard  approach 
(i.e.,  the  original  application  of  equation  35)  would  yield  the  solution  valid 
to  the  second  order,  the  next  iteration  would  yield  the  solution  valid  to  the 
third  order,  and  so  on.  To  summarize,  the  order  of  errors  in  the  repeated 

applications  of  (36)  would  be  4 ,  6 ,  8 . while  this  order  in  the  standard 

iterative  process  symbolized  by  (35)  would  be  3,  4,  5 . 

From  the  theoretical  standpoint,  the  most  revealing  Is  the  third  case.  As 
we  have  seen,  the  original  formulation  (23a)  is  valid  to  within  the  third  order. 

P 

However,  if  <5x  contains  finite  quantities,  the  standard  form  of  normal 

equations  (35)  is  virtually  useless;  It  is  valid  only  to  the  first  order,  since 

the  first  neglected  term  in  (23a).  namely  -(<5x  ) nf  AlAu®),  is  of  the  second 

s  pa 

order,  the  same  as  the  entire  left-hand  side  of  (35).  Thus,  not  only  does  (35) 
entail  errors  two  orders  of  magnitude  greater  than  errors  in  (23a),  but,  in 

p 

general.  It  does  not  converge.  (If  Ax  should  belong  to  an  "In-between" 
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category  discussed  in  Section  3.2,  these  statements  would  be  comaensurabiy 
softened.)  Such  facts  are  at  the  root  of  the  standard  assumption  that  not  only 
and  thereby  (du**),  but  also  dx^  and  thereby  should  be  small.  We 

conclude  by  stating  that  the  standard  approach  cannot  accommodate  the  third 
case,  while  the  new  formula  (23a)  accommodates  it  to  the  third  order,  and  the 
simplified  formula  (36)  accommodates  it  to  the  second  order. 

As  a  matter  of  interest,  we  examine  the  possibility  of  d(du*)  being  larger 
than  second-order.  Out  of  necessity,  we  assume  that  d(du**)  contains  quantities 
whose  magnitude  is  somewhere  between  the  first  and  the  second  orders,  in  the 
sense  that  o[d(Au**)  A{du^)]  is  of  the  third  order  (instead  of  the  fourth  order 
as  stipulated  earlier).  Clearly,  the  neglected  quantities  o[A(Au  )  A(Au^)] 
cannot  be  of  the  second  order  since  this  would  mal<e  even  the  full  equation  {23a) 
valid  only  to  the  first  order.  The  dots  in  equations  (17b)  and  (19b)  now 
represent  third-  (instead  of  fourth-)  and  higher-order  quantities.  The  same  can 
be  said  about  the  terms  a,  b,  and  c,  and  thereby  about  (22)  and  about  errors  in 
(23a).  Due  to  the  modification  of  the  assumption  regarding  A{Au*),  the  accuracy 
of  the  three  formulas  under  consideration  (equations  35,  36,  and  23a)  is 
somewhat  compromised  in  the  first  and  the  second  cases,  but  the  effect  is  not 
qualitatively  profound.  However,  in  the  third  case  the  standard  formula  (35)  is 
virtually  useless.  The  simplified  formult,  (36).  although  compromised  to  a 
certain  extent,  leads  to  useful  results,  while  the  full  formula  (23a)  is  as 
effective  as  the  standard  formula  (35)  would  be  under  the  standard  assumptions. 

r 

(Here  again,  if  dx  should  belong  to  an  "in-between"  category,  all  three 
formulas  would  give  commensurably  better  results.) 


We  note  that 
evaluated  at  ( P) . 


(23a)  could  be  reformulated  with  the  partial 
Since  at  (P)  we  have 


i3a 


Vr  ^ 


derivatives 


{37a) 


In  l<eeplng  with  the  stipulated  error  level  we  can  write  (23a)  as 


[(a^^)  -  (<5Xg)  (n|^)I  A{Au®)  =  (A®)  (dx^)  .  (37b) 

However,  the  original  equation  (23a)  is  preferable  to  (37b)  since  It  displays 

s  s 

the  fact  that  the  sets  D.  and  <P,  can  be  stored  and  treated  as  constant, 

pa  par 

and  implicitly  suggests  In  which  manner  and  under  what  conditions  It  can  be 
simplified.  We  have  seen  such  simplifications  In  (35)  and  (36).  Although  the 
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full  formula  (23a)  is  appealing  from  the  theoretical  standpoint,  computer  run¬ 
time  and  storage  considerations  favor  the  simplified  formula  (36). 


After  obtaining  the  final  parametric  values  u**-u*+Au*,  one  can  compute 

r  -r  ^ 

the  desired  quantities,  such  as  Ax  and  thereby  dx  ,  as  well  as  the  tensor 

needed  to  form 
a 


a 


^|3  ®sr 


(38) 


which,  in  turn,  can  be  used  In  the  variance-covariance  propagation  and  In  the 
weight  propagation.  In  paralleling  the  discussion  of  Section  3.1  pertaining  to 
a  linear  model,  we  state  that  the  tensor  a  obtained  from  a.  via  the 
tensor  equation  a  ^a^^=d^,  corresponds  to  the  variance-covariance  matrix 
of  adjusted  parameters,  while  the  tensor  a^^  Itself  corresponds  to  the  weight 
matrix  of  these  parameters.  Further,  the  variance-covariance  and  the  weight 
matrices  of  adjusted  observations  correspond  respectively  to 


-,rs  -rr  -aB  ::3 
g-  =  A^  a  A^  , 


si  ^  ^jr 


(39a, b) 


As  has  been  mentioned  at  the  close  of  Chapter  2,  the  above  tensors  describe 
"tangential",  or  first-order,  properties  of  least  squares  estimates.  In 
reference  to  [Blaha,  1984],  the  variance-covariance  and  the  weight  matrices  of 
residuals  correspond  respectively  to  the  tensors 


-„rs  rs 
g  =  g 


.,rs 


sr 


■  S. 


sr 


(40a, b) 


4 . 2  Additional  Notes 

We  begin  this  discussion  by  demonstrating  that  an  eventual  re-appi Ication 
of  (23a)  is  straightforward  and  without  any  pitfalls.  We  observe  that  such  a 
re-application  would  be  equivalent  to  using  (23a)  in  conjunction  with  a  new 
point  (P)  situated  closer  to  P  than  the  original  point  (P)  obtained  via  the 

V 

Initial  iteration.  The  quantity  (Au  )  would  then  represent  the  model-surface 

coordinate  differences  between  the  new  point  (P)  and  the  Initial  point  P;  (<5x  ) 

would  similarly  refer  to  the  new  point  (P).  This  situation  would  change  nothing 

in  the  development  nor  in  the  results  represented  by  equation  (23a)  Itself.  In 

particular,  the  utilization  of  (Au„)  from  (8),  which  in  Fig.  1  gave  rise  to  the 

P 

model-plane  point  P^  and  thereby  to  the  model-surface  point  (P).  was  essentially 
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a  matter  of  convenience.  If  3uch  a  model-plane  point  were  found  closer  to 
(produced  e.g.  by  the  Initial  iteration  of  Section  3.2),  It  would  represent  a 
welcome  shortcut,  but  none  of  the  pertinent  relationships  from  (12a)  on  would 
change.  And  If  It  could  be  found  at  the  exact  location  of  P^.  the  right-hand 
side  of  (23a)  would  be  zero  by  virtue  of  (7),  which  would  yield  A(Au*)-0.  The 
same  considerations  are  valid  also  with  regard  to  the  simplified  formula  (36). 

As  we  have  seen,  the  latter  differs  from  the  full  formula  (23a)  by  the  error 
level  due  to  the  truncated  terms,  but  only  In  the  case  of  large  residuals. 

We  now  focus  our  attention  on  the  second  case  (with  small  residuals) 

I* 

represented  by  (36).  and  examine  whether  advantage  can  be  taken  of  the  set 

in  forming  (A^)  and  (<5x^).  as  opposed  to  computing  these  tensors  from  the 

nonlinear  model.  If  we  Include  only  the  term  containing  0^^  In  (16)  and 

neglect  the  subsequent  terms,  the  error  In  (A^)  will  be  of  the  second  order, 

as  will  be  the  error  In  (a.  )  formed  through  (23b).  Accordingly,  the  error  on 

pa 

the  left-hand  side  of  (36)  will  be  of  the  fourth  order  and  negligible. 

Similarly,  If  we  Include  only  the  term  containing  o'*  In  (18)  and  neglect  the 

r  ^  r 

subsequent  terms,  the  error  in  (Ax  )  and  thereby  also  in  (<5x  )  given  by  (20b) 
will  be  of  the  third  order,  resulting  only  In  a  fifth-order  error  on  the  left- 
hand  side  of  (36) . 

However,  the  right-hand  side  of  (36)  offers  a  different  picture.  In 
forming  (Aa)(i5x  ),  where  (Af)  contains  second-order  errors  and  the  flrst- 

p  3  p 

order  quantity  ((Sx^)  contains  third-order  errors,  we  observe  that  either  kind  of 

errors  introduces  third-order  errors  into  the  product.  Thus,  to  keep  with  the 

r  r 

fourth-order  error  level  implied  by  (36),  both  (A  )  from  (16)  and  (Ax  )  from 

r 

(18)  would  have  to  include  the  terms  with  We  note  that  the  order  of 

magnitude  of  the  remaining  truncation  errors  In  these  two  quantities  is 
independent  of  the  closeness  of  (Au**)  to  Au®,  which  would  preclude  a  re- 
application  of  (36). 

In  the  case  of  a  single  application  of  (36),  only  (16)  would  entail 

significant  computing  effort  since  (18)  could  take  advantage  of  the  intermediate 

results  fi'^-(Au^)  and  (Au^)(Au^).  If  one  chose  to  pursue  this  avenue,  it 
ap  g  _ 

would  then  be  quite  easy  to  adopt  the  full  formula  (23a)  where  ®^^^(Au  )  would 
be  known  from  an  early  intermediate  result.  However,  passing  from  (36)  to  (23a) 
would  not  decrease  the  overall  error  level  (fourth-order).  In  realizing  that 
<6^  is  not  needed  In  (38),  we  might  wish  to  avoid  computing  and  storing  these 
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third-order  partial  derivatives  altogether.  Accordingly,  (A^)  as  well  as 
(Ax*’)  and  thereby  (dx*')=fix*'-(Ax*')sXQ-(x*') ,  where  (x*’)-x'^(  (u*) ) ,  would  be 
cowputed  from  the  nonlinear  model  as  In  the  standard  adjustment  approach.  This 
choice  is  further  supported  by  a  potential  need  to  re-apply  the  algorithm  (36). 


It  is  instructive  to  compare  the  development  presented  in  this  study  with 

that  of  [Blaha,  1987].  Although  both  cases  rely  on  the  precise  relation  (7). 

the  latter  resulted  In  a  less  efficient  algorithm  due  to  the  simplification 

**s 

which  we  shall  now  describe.  The  tensor  A„  was  utilized  In  the  form 

P 

A^+AAq  as  In  (13a).  On  the  other  hand,  the  vector  Sx  was  substituted  for 
P  P  s 

by  the  right-hand  side  of  (14),  but  only  for  the  contraction  with  A^g^^.  not 

for  the  contraction  with  AA® .  Accordingly,  the  relation  (7)  became 


A^  g  dx*'  -  A®  g  Ax*^  AA®  g  dx*^ 

(9  ^sr  $  ^sr  0  ^sr 


0  . 


The  first  term  above  was  (Au„)  by  virtue  of  (8),  and  the  second  term  was 

p 

-a^  Au^-A^g  v'*.  as  Is  confirmed  through  (12a, b).  The  resulting  equation  was 
0r  $  sr 

contracted  with  a  yielding 


A  «  /A  ot\  .3  r  a$  ..s  .-r 

Au  =  (Au  )  -  a  ^  A^  gg^  V  a  AA^  g^^  dx  , 

which,  except  for  minor  notational  differences,  was  presented  in  [Blaha,  1987] 
as  equation  (28). 


The  solution  in  this  reference  proceeded  as  If  Au**  were  split  into 
( Au** )  ♦^ A ( Au® )  seen  in  (15),  but  only  for  the  left-hand  side  above,  not  for  the 
formation  of  v*^  by  (12c),  nor  for  the  formation  of  dx*^  by  (14)  with  Ax*^  given 
by  (12a-c).  The  values  Au®  implied  In  v*'  and  Ax*'  were  approximated  by  (Au®)  , 
i.e.,  the  terms  beyond  the  first  on  the  right-hand  side  of  the  above  equation 
were  treated  as  if  A(Au®)  were  zero.  This  introduced  errors  in  an  Initially 
exact  equation,  which  were  an  order  of  magnitude  greater  than  the  terms 
neglected  In  the  present  study.  For  example.  In  the  situation  called  here  the 
second  case  the  thus  Introduced  errors  would  be  of  the  third  order  Instead  of 
the  current  fourth  order.  Such  simplifications  downgraded  the  convergence 
characteristics  of  [Blaha,  1987]  essentially  to  the  level  of  the  standard 
approach . 


The  results  as  well  as  intermediate  formulas  in  this  study  could  be 
prece.nt<»d  in  the  familiar  matrix  notation.  In  the  case  of  simply-  and  doubly- 
indexed  quantities,  matrix  transcriptions  would  present  no  difficulty  since 
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tensor  contractions  and  matrix  multiplications  amount  to  the  same  operations, 
provided  the  matrices  are  arranged  in  the  proper  sequence  and  form  (direct  or 
transposed).  With  regard  to  three-  and  higher-dimensional  arrays,  matrix 
manipulations  could  be  extended  upon  imagining  that  ordinary  matrices  are 
stacked  one  behind  the  other,  that  groups  of  matrices  are  stacked  one  behind  the 
other,  etc.  However,  we  prefer  working  with  the  original  Indexed  quantities, 
which  allows  more  flexibility  in  derivations  as  well  as  computations,  and  which, 
in  any  case,  is  at  the  root  of  the  matrix  operations  themselves  (restricted  to 
one  or  two  indices).  For  the  sake  of  Interest,  Appendix  D  translates  the  main 
outcome  of  this  study  into  matrix  notation. 


5.  SOIMARY  AND  CONCLUSIONS 


A  convenient  approach  for  resolving  nonlinear  least-squares  probleas 

consists  In  using  an  Isomorphic  geometrical  setup  with  tensor  structure  and 

notation.  Such  a  link  Is  highlighted  by  the  consideration  that  If,  In  the  n- 

dlmenslonal  observational  space,  the  metric  tensor  g  corresponds  to  the  weight 

-r 

matrix  P,  and  if  the  contravarlant  components  dx  of  a  displacement  symbolized 
by  Sx  correspond  to  a  set  of  estimated  observational  errors,  l.e.,  to  a  set 
of  negative  residuals  -V.  then  the  least-squares  criterion 

P  V  =  minimum 


corresponds  to  the  quadratic  form 

—  3  —  r  ~  2 

<5x  g  6x  s  Ss  =  minimum  , 
sr 

where  <53  is  the  magnitude  of  the  displacement.  Should  5x  identify  an 
unknown  displacement  from  the  u-diraenslonal  model  surface  to  the  point  Q  given 
in  the  surrounding  observational  space  through  the  actual  observations,  the 
minimum-distance  property  restricts  <5x  to  the  segment  PQ  of  a  straight  line 
dropped  orthogonally  from  Q  onto  the  model  surface,  where  P  marks  the 
intersection.  This  situation  is  illustrated  schematically  in  Pig.  1.  Since 
orthogonal  relations  are  readily  exploited  In  a  geometrical  environment,  the 
least-squares  method  --  unlike  any  other  method  --  is  particularly  suited  for 
geometrical  analogy  and  analysis. 


Among  the  basic  correspondences  between  adjustments  and  geometry,  the 

number  of  observations,  n,  and  the  number  of  parameters,  u,  correspond 

respectively  to  the  dimensionality  of  the  observational  space  and  of  the  model 

surface.  The  variance-covariance  matrix  of  observations,  E,  corresponds  to  the 

rs 

associated  metric  tensor  g  ,  while  the  weight  matrix  of  observations,  adopted 

as  P»E  corresponds  to  the  metric  tensor  g  Since  I  and  thus  also  P  are 

sr 

constant,  the  observational  space  is  endowed  with  a  coordinate  system  {x  }  such 
that 


g  =  constant 
sr 


rs 

g 


constant 


b  r 

The  set  L  of  actual  observations  corresponds  to  the  set  x^ 

space  coordinates  describing  the  point  Q.  All  possible  sets 


of  observat i onai - 
L*  of  adjusted 
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observations  (subject  to  no  criterion)  correspond  to  the  Gauss  form  of  the  model 
surface  endowed  with  a  coordinate  system  (u®}; 

r  r ,  a.  ,  - 

X  =»  X  (u  )  .  r  =  1.2 . n  .  a  =  1,2 . u  . 


The  final  set  of  adjusted  parameters,  X  ,  corresponds  to  a  particular  set 

u®  of  model-surface  coordinates  describing  the  least-squares  point  P.  The  set 

of  initial  parameters,  X°.  corresponds  to  the  set  u®  of  model-surface 

coordinates  describing  the  Initial  point  P.  The  final  set  of  parametric 

corrections,  X,  then  corresponds  to  Au®=u®-u®;  these  quantities  are  assumed 

to  be  small  (termed  first-order).  The  final  set  of  adjusted  observations,  L  , 

r  r  ot 

corresponds  to  a  particular  set  x  =x  (u  )  of  observational -space  coordinates 

describing  the  least-squares  point  P.  Similarly,  the  set  P(X®)  corresponds  to 

the  set  x''=x*'(u®)  of  observational-space  coordinates  describing  the  initial 

point  P.  The  set  of  negative  constant  terms.  -L  =  l“-P(x”),  corresponds  to  the 

r  r  p 

contravarlant  vector  <5x  ,  while  the  set  of  negative  residuals, 

^  Q  o 

Da  ”  r  r  r 

-V  =  L  -P(X  ).  corresponds  to  the  contravarlant  vector  «5x  ■'X^-x  .  The 

Initial  design  matrix  A,  which  in  standard  observation  equations  V=AX+-L  relates 

the  parametric  corrections  to  the  residuals,  corresponds  to  the  design  tensor 

A^=3x'^/au®  evaluated  at  P.  On  the  other  hand,  the  standard  adjustment 

approach  does  not  have  equivalents  of  Q**.  and  ,  which  form  three-  and 

ap  apy 

four-dimensional  arrays,  respectively,  and  contain  the  second-  and  the  third- 

r  Cl 

order  partial  derivatives  of  x  with  respect  to  u  ,  evaluated  at  P. 

The  geometrical  approach  is  based  on  a  direct  exploitation  of  the  relation 


A®  g^^  <5x^  >  0  ,  (41) 

p  s  r 

where  A®  represents  the  design  tensor  evaluated  at  the  least-squares  point 
_  P 

P,  and  equation  (41)  Itself  represents  the  orthogonality  condition  at  P. 

The  outcome  of  the  geometrical  development  is  considered  in  two  methods,  called 
geometrical  and  extended  geometrical.  It  is  contrasted  to  the  standard  method 
treating  nonlinear  models  In  a  linearized  form.  The  algorithms  associated  with 
all  three  methods  are  presented  below  In  the  form  of  the  first  Iteration,  and  in 
t.ie  form  of  the  second  and  subsequent  Iterations. 


All  the  formulas  are  written  in  tensor  notation.  In  the  case  of  slmply- 
and  doubly-indexed  quantities,  matrix  transcriptions  would  present  no  difficulty 
since  tensor  contractions  and  matrix  multiplications  amount  to  the  same 


operations,  provided  the  matrices  are  arranged  in  the  proper  sequence  and  form 
(direct  or  transposed).  With  regard  to  three-  and  higher-dimensional  arrays, 
matrix  manipulations  could  be  extended  upon  imagining  that  ordinary  matrices  are 
stacked  one  behind  the  other,  that  groups  of  matrices  are  stacked  one  behind  the 
other,  etc.  However,  working  with  the  original  Indexed  quantities  allows  more 
flexibility  in  derivations  as  well  as  in  computations. 


In  tensor  notation,  the  Initial  matrix  of  normal  equations  corresponds  to 


the  model-surface  metric  tensor  a^^  at  the  initial  point  P,  and  the  initial 
right-hand  side  of  normal  equations  corresponds  to  the  model-surface  covariant 


vector  A  5x  at  P,  where 
p  s 


Aa  g  A 
$  sr  a 


/9a 


dx 


g-_  <5x 
sr 


The  parametric  corrections  stemming  from  the  first  iteration  are  symbolized  by 
(Au**),  and  they  give  rise  to  an  updated  point  (P).  The  latter  is  described  by 
the  model-surface  coordinates  (u**)  ”U^+ (Au**) .  The  quantities  belonging  to  (P) 
will  likewise  be  written  in  parentheses.  The  parametric  corrections  obtained  in 
the  second  iteration  will  be  denoted  A(Au*),  and  they  will  give  rise  to  a  new 
updated  point  determined  via  (u®)>A(Au“).  The  notation  used  in  conjunction  with 
the  second  iteration  will  be  retained  also  for  any  further  iterations. 

cc  r 

standard  method.  Under  the  assumption  that  both  sets  Au  and  <5x  contain 
small  quantities  (first-order),  the  first  iteration  in  the  standard  method  reads 


reoresenting  the  initial  normal  equations.  The  second  and  further  iterations 
follow  the  same  principle: 


(a^^)  A(Au®)  =  (A®)  (<5Xg)  .  (42b) 

As  is  illustrated  by  the  example  in  Appendix  A,  the  convergence  properties  of 
this  method  may  be  significantly  impaired  by  an  increased  size  of  observational 
errors . 


Geometrical  method.  Under  the  same  assumption  as  above  (both  sets  Au**  and 
<5x^  contain  small  quantities),  the  first  iteration  utilizes  the  same  formula  as 
its  standard  counterpart,  namely 
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However,  the  second  and  further  Iterations  proceed  according  to 


[(a^^)  -  (fiXg)  A(Au®)  =  (A®)  (5Xg)  .  (43b) 

representing  the  modified  normal  equations.  The  triply-indexed  quantity  0^^, 
formed  by  second-order  partial  derivatives  of  the  observables  with  respect  to 
the  parameters,  is  evaluated  only  at  the  Initial  point  P.  The  formulas  (43a, b) 
are  relatively  simple  yet  efficient,  and  are  suitable  for  practical  applications 
where  the  convergence  of  the  standard  method  Is  problematic. 


d 

Extended  geometrical  method.  Although  the  assumption  regarding  Au  is 
unchanged,  this  method  Is  tailored  for  dx^  containing  relatively  large 
quantities,  for  which  the  first  Iteration  reads 


Compared  to  the  geometrical  method,  the  current  algorithm  Is  seen  to  utilize 
second-order  partial  derivatives  and  to  give  rise  to  the  modified  normal 
equations  already  In  its  first  Iteration.  The  formula  for  the  second  and 
further  iterations  Is  given  as 


representing  the  modified  normal  equations  at  updated  stages.  Here  use  Is  made 

of  containing  third-order  partial  derivatives.  This  quantity  Is 

iSar  3 

evaluated  only  at  P,  similar  In  this  respect  to  0“  .  We  note  that  the 

^  s 

quantity  inside  the  brackets  of  (44b)  could  be  replaced  by  ~  (<5x^ )  (Q^^) , 


s 


/3a' 


where  (Q.  )  would  represent  the  second-order  partial  derivatives  evaluated  at 
P« 

an  updated  point. 


After  obtaining  the  final  set  of  parametric  values  u®-u**<-Au®,  one  can 

— r  r  ~r 

compute  the  desired  quantities  A  ,  Ax  ,  and  Ax  either  from  the  nonlinear 

— r 

model  or  by  using  the  Taylor  series  expansion.  The  tensor  A^  Is  needed  to 
form 


®3r 


-aB 

which  corresponds  to  the  weight  matrix  of  adjusted  parameters.  The  tensor  a 

obtained  from  a.  through  the  relation  a**^a„  =*6*  then  corresponds  to  the 
pa  P7  T 

variance-covariance  matrix  of  adjusted  parameters.  Further  we  have 


-,rs  -r  -a0  -s 
g  =*  A  a  Ao  . 
^  a  0 


Kt  *  hi  'j.- 
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which  correspond  respectively  to  the  variance-covariance  and  che  weight  Matrices 
of  adjusted  observations.  The  above  tensors  are  expressible  by  aeans  of 
orthonoraal  vectors  tangent  to  the  aodel  surface  at  P.  Accordingly,  they 
describe  "tangential",  or  first  order,  properties  of  least-squares  estlaates. 
Finally,  the  tensors 


i^rs 


rs 


-  ,ra 
g' 


’sr 


sr 


-  g* 

*sr 


correspond  respectively  to  the  variance-covariance  and  the  weight  Matrices  of 
residuals . 


A  significant  contribution  to  the  understanding  of  the  least-squares 
adjustment  theory  through  differential  geometry  has  been  made  by  Teunlssen 
[1985].  The  geometry  of  a  nonlinear  adjustment  is  treated  In  his  Chapter  IV, 
the  aaln  topic  of  which  Is  Gauss'  Iteration  method.  This  method  characterizes 
the  standard  adjustment  approach  as  Is  conflraed,  for  example,  by  the  text  on 
page  109  stating  that  "at  each  Iteration  step  the  observation  point  is 
orthogonally  projected  onto  a  new  tangent  space".  Such  an  algorithm  is 
equivalent  to  the  solution  of  the  type  (42a, b)  above,  which,  in  terms  of  Pig.  1, 
would  result  In  the  projection  of  the  point  Q  onto  the  aodel  plane  passing 
through  the  Initial  point  P,  followed  by  the  projection  of  Q  onto  a  new  aodel 
plane  passing  through  an  updated  point  (P),  etc. 

In  the  sane  vein,  page  109  (Ibid.)  states  that  Gauss'  iteration  method 
"consists  of  successively  solving  a  linear  least-distance  adjustment  until  the 
necessary  condition  of  orthogonality  is  fulfilled".  This  condition,  equivalent 
to  our  equation  (41),  Is  fulfilled  in  the  standard  method  as  a  by-product  of  the 
projections  discussed  above.  By  contrast,  the  geometrical  approach  actively 
seeks  to  fulfill  It  at  every  step.  A  one-step  solution  producing  the  least- 
squares  point  P  directly  Is  hindered  only  by  the  necessity  to  truncate  certain 
terns,  but  not  to  the  extent  of  making  the  entire  model  linear  (see  the  above 
equations  43b  and  44a, b).  The  matrix  of  modified  normal  equations  generated  In 
the  process  Is  positive-definite,  similar  in  this  respect  to  the  matrix  of 
normal  equations  In  the  standard  method. 

Although  the  behavior  of  the  standard  adjustment  algorithm  has  been 
analyzed  and  described  In  detail  by  Teunlssen  [1985],  the  behavior  of  the 
geometrical  algorithm  developed  herein  allows  ample  opportunity  for  further 
analysis.  A  study  comparing  the  geometrical  approach  with  alternative  methods 
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would  be  worthwhile  for  the  sake  of  theory  as  well  as  for  practical  reasons 
(conputatlonal  burden,  rate  of  convergence,  etc.).  Another  worthwhile  effort 
would  be  to  study  the  effect  of  the  model's  nonlinearity  on  the  statistical 
properties  of  least-squares  estimators,  as  is  outlined  In  Section  6  of  Chapter 
rv,  ibid.  The  need  for  further  research  is  also  underscored  by  the  fact  that 
the  above  discussion  has  touched  merely  on  one  class  of  adjustment  problems,  the 
parametric  adjustment.  Nonlinear  versions  of  the  general  method  (with  or 
without  constraints),  rank-deficient  systems,  etc.,  all  await  a  systematic 
geometrical  treatment. 

Indeed,  as  the  Epilogue  to  his  Chapter  IV.  Teunlssen  [1985]  indicates  that 
scientists  are  only  on  the  brink  of  understanding  the  complexity  of  nonlinear 
adjustment.  He  further  states:  "Unfortunately,  one  can  seldom  extend  the 
elegant  formulation  and  solution  techniques  from  linear  to  non-Jlnear 
situations.  For  most  non-linear  problems  one  will  therefore  have  to  recourse. 

In  practice,  to  methods  which  are  iterative  In  nature”.  In  this  light,  the 
present  study  constitutes  a  contribution  in  accelerating,  although  not 
completely  eliminating,  the  iterative  process. 

Encouraging  results  have  been  obtained  in  the  numerical  example  presented 
in  Appendix  A.  illustrating  convergence  properties  of  a  third-order  polynomial 
in  four  variables.  Although  the  standard  method  converged  slowly  In  one  of  four 
analyzed  cases  and  diverged  in  two  others,  the  geometrical  method  converged  in 
two  and  three  iterations,  respectively.  The  extended  geometrical  method  further 
reduced  the  number  of  iterations  from  chree  to  two.  I u  lo  uxpected  that  in  most 
nonlinear  cases  the  presence  of  second-order  partial  derivatives  will  translate 
into  two  iterations  in  the  geometrical  method  as  compared  to  several  Iterations 
needed  by  the  standard  method. 

In  spite  of  the  accelerated  convergence  of  the  geometrical  approach,  the 
fact  remains  that  the  methods  discussed  herein,  as  well  as  other  methods  used  in 
physical  sciences  for  adjusting  nonlinear  models,  are  Iterative  in  nature. 
However,  one  can  conceive  of  a  geometrical  scheme  where  the  determination  of  the 
least-squares  point  P  can  be  made  arbitrarily  accurate  in  one  single  step  upon 
increasing  the  number  of  meaningful  terms  In  the  Taylor  series  expansions.  Such 
an  approach  would  have  to  be  fundamentally  different  from  the  one  described  in 
this  study,  where  an  Increase  in  the  number  of  meaningful  terms  has  certain 
limitations.  In  particular,  if  the  observational  errors  are  relatively  small. 
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the  current  geometrical  approach  cannot  take  advantage  of  the  third-  and  higher- 
order  partial  derivatives  of  the  observables  with  respect  to  the  parameters, 
whose  contribution  would  be  Invalidated  by  the  truncation  errors  of  comparable 
magnitudes.  Thus,  developing  a  theoretically  distinct  geometrical  setup  along 
these  lines  Is  yet  another  challenge  In  the  field  of  the  nonlinear  least-squares 
adjustment. 


34 


APPENDIX  A 


NUNERICAL  EXAMPLE 


The  nonlinear  parametric  adjustment  model  is  represented  by  six 
observables,  each  expressed  by  a  third-order  polynomial  in  four  variables. 

Thus,  the  number  of  observations  is  six,  n=*6,  and  the  number  of  parameters  is 
four,  u=4.  A  simpler  adjustment  problem  featuring  a  second-order  polynomial  was 
presented  by  Bessette  [1987],  who  utilized  the  algorithm  described  in  [Blaha, 
1987],  Although  that  algorithm  lowered  the  number  of  scalar  multiplications  in 
each  iteration  as  compared  to  the  standard  method,  it  did  not  accelerate  the 
rate  of  convergence  (see  Section  4.2).  It  has  been  superseded  by  the  algorithm 
developed  here,  which  has  been  programmed  and  verified  independently  on  two 
different  computers. 


In  terms  of  the  four  parameters  u  ,  a=l,2,3,4,  a  general  observable  denoted 
X  is  expressed  by 


Ot  (X  S  CL  S  T 

'  "as  “  “  '  •a0r  "  “  "  • 


where  a  ranges  from  1  to  4.  (9  ranges  from  a  to  4,  and  r  ranges  from  0  to  4.  The 
numerical  values  of  the  four  c-coef f iclents  (linear),  the  ten  k-coefflclents 
(quadratic),  and  the  twenty  m-coeff Iclents  (cubic)  for  all  six  observables  are: 
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-2.0, 

*1.7,  *0.8, 

+  1.0, 

-2.1. 

*2.2, 

*1.2. 

-1.9, 

+  1.7, 

-1.4. 

-0.8, 

*0.8, 

*1.0. 

-1  - 1  , 

-0,9,  +1.8, 

-1.4, 

-0.9, 

*2.1  . 

*1.6, 

*1.0, 

-2.1. 

*1.1: 

-1.7, 

*2.2, 

*1.6, 

-1.6, 

*1.6, 

-1.2, 

*2.1, 

*2.3, 

-2.1, 

*1.2; 

-2.2, 

-0.8, 

+1.0,  -1.2. 

+0.9, 

*1.3, 

-0.8, 

+2.1, 

*1.9, 

*1.2, 

+0.9. 

-1.5, 

-16, 

-0.8. 

-1.0, 

+0.9,  -0.9, 

+  2.1, 

*0.8, 

-1.3. 
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+0.9, 

+2.1, 

-2.1,  +1.8: 

+  1.2, 

+  1.7, 

+2.2, 

-1.5, 

+  1.9, 

+  1.3, 

-2.2, 

+  1.7, 

+  1.0, 

+  2.3; 

+  1.1, 

-1.2,  -0.9, 

+  1.3, 

+  1.5, 

+  1.0, 

-0.8, 

-2.2, 

+0.8, 

+  1.7, 

-1.2, 

-0.3, 

+  1.1, 

-2.0, 

+1.7,  -0.9, 

+0.8. 

+1.3, 

-1.8, 

-1.1. 

+  1.7, 

+  1.2, 

+1.7,  -1.6; 

+2.1, 

+  1.6, 

-1.1, 

+  1.9, 

+  1.6, 

+  1.7, 

+2.0, 

+  1.2, 

+  1.6, 

-1.7; 

+  1.3, 

+2.1,  -1.0, 

-0.8, 

+1.2, 

-2.1, 

+2.2, 

-1.6, 

-1.8, 

-0.8, 

+  1.9, 

+  1.2, 

-0.9, 

+  1.2, 

+2.0,  -1.7, 

+1.1, 

+1.0, 

+0.8, 

-2.2. 

The  positive-definite  variance-covariance  aatrlx  of  observations 

PS 

corresponds  to  the  associated  aetrlc  tensor  g  ,  chosen  as 


^1.3266  +0.5254 

+0.3669 

+0.1572 

-0.2667 

-0.0549 

+1.3457 

-0.0281 

+0.5924 

-0.1448 

+0.5521 

+1.4143 

+0.5144 

+0.3600 

-0.5204 

. 

+1.4629 

-0.5604 

+0.4877 

. 

. 

. 

+1.3930 

-0.1034 

+1.3500 

where  the  dots  denote  the  symmetric  elements.  The  positive-definite  weight 

matrix  of  observations  corresponds  to  the  metric  tensor  g  computed  from 
.  .  sr 

g  Initial  set  of  parameters  Is  represented  by 

u®  >  +1.399,  -1.201,  +1.299,  +1.099  , 

describing  the  Initial  point  P  by  model-surface  coordinates.  The  observational- 

space  coordinates  x  of  this  point  are  determined  from  the  model,  as  are  the 

®  r  r  r 

sets  of  partial  derivatives  A  ,  0  and  ®  „  (the  last  set  is  constant 

<x  otp  ^pT  P 

here,  and  any  further  set  Is  zero).  Finally,  the  elements  of  x  representing 
error-free  observations  are 

X**  -  +4.80,  -4.79,  +33.44,  +0.33,  +22.47,  -11.69  . 

These  values  correspond  to  the  following  set  of  final  parameters: 

u®  >■  +1.4,  -1.2,  +1.3,  +1.1  . 

The  example  has  been  generated  In  four  versions  which  differ  from  one 
another  by  the  size  of  the  final  residuals.  Each  version  Is  treated  by  three 
methods  summarized  in  the  conclusion:  standard,  geometrical,  and  extended 
geometrical.  A  basic  version  has  been  generated  with  excessively  large  final 
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—r 

residuals  grouped  in  the  tensor  -i5x  ; 

<Sx^  «  +3.22,  +22.55,  -10.59,  -1.28,  +24.27,  +14.68  . 

-r 

This  version  gives  rise  to  a  total  of  four,  upon  scaling  the  values  In  <5x  by 
a  factor  c,  where 

c  =  0.001:  0.01;  0.1:  1. 

The  residuals  In  the  first  version  (c=0.001)  can  be  considered  first-order,  or 
essentially  differential,  as  required  by  the  standard  methodology.  On  the  other 
side  of  the  spectrum,  the  residuals  in  the  fourth  version  (c=l)  are  so  large 
that  this  case  is  presented  mainly  for  the  sake  of  theoretical  interest. 

The  iterative  process  is  terminated  in  each  method  and  each  version  when 

ot  “8 

errors  In  all  of  the  parameters  u  ,  a*l,2,3,4,  are  less  or  equal  to  5x10 

(K  2 

Such  errors  remaining  in  u  are  listed  below  only  for  u  ,  for  which  they  are  in 

4 

general  the  worst.  In  the  standard  method,  the  errors  In  u  are  very  similar  to 
2 

those  in  u  (Including  the  sign),  but  their  magnitude  Is  slightly  larger.  A 
statement  of  this  kind  applies  also  to  the  geometrical  method  for  small  values 
of  c  (represented  by  c=0.001  and  c=0.01).  The  reverse  Is  generally  true  for 
large  values  of  c  (represented  by  c=0.1  and  c»l)  in  the  geometrical  method,  and 
for  all  four  c's  In  the  extended  geometrical  method. 

Results  for  u  In  the  standard  method.  In  the  geometrical  method,  and  In 
the  extended  geometrical  method  are; 


c  =  0.001  . 

.  stand. ; 

-4x10  ®, 

+3xl0"®  ; 

geom . ; 

-4xl0”®, 

+2xl0”^^ 

c  =•  0.01 

.  stand . : 

-3xl0'^, 

+2xl0'®. 

-3x10  . 

+3xl0"®  : 

geom . : 

-3xl0'®, 

+2xl0”^  . 

c  =  0. 1 

.  stand. ; 

-4 

-2x10  . 

+2x10”'*, 

\ 

X 

o 

1 

'4 

^3x10  .  diverges 

geom. : 

-2xl0"*, 

+lxl0”^. 

-IxlO 

; 

ext .  geom. : 

-3xl0"^, 

-Ixl0"^° 
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stand. : 

-2xl0"®, 

+2x10”^, 

-1  0 

-2x10  ,  +1x10  .  diverges 

geom .  : 

-2xl0”^. 

-ixio"®. 

+lxlO~®  ; 

geom . : 

-Ixl0“®, 

-2xl0"^° 

The  above  results  and  further  coaputer  runs  warrant  the  following  coaaents: 

2 

1)  Results  substantially  worse  than  those  listed  for  u  have  been  noticed 

4  -11 

for  C‘0.001  In  the  case  of  the  parameter  u  In  the  geoaetrlcal  aethod  (+1x10 

-12 

as  opposed  to  +2xio  ).  However,  such  exceptions  do  not  detract  from  the 
overall  pattern  seen  above. 

2)  The  passage  from  the  geometrical  aethod  to  the  extended  geometrical 
method,  l.e.,  the  inclusion  of  second-order  partial  derivatives  in  the  first 
iteration  and  of  third-order  partial  derivatives  In  the  second  Iteration,  has  no 
appreciable  effect  on  the  adjustment  with  small  residuals  (c=0.001  and  c=0.01). 
On  the  other  hand,  in  case  of  large  residuals  (c=0,l  and  c=l)  the  effect  Is  seen 
to  be  great  and  the  number  of  Iterations  is  reduced  from  three  to  two. 

3)  If,  in  the  extended  geometrical  method  for  c=»0.1  and  c=l ,  the  first 
Iteration  is  adopted  from  the  standard  method  (l.e.,  if  the  second-order  partial 
derivatives  are  excluded  from  that  iteration),  the  total  number  of  iterations 
will  grow  from  two  to  three.  As  an  example  in  the  case  c*l ,  the  errors  in  the 
third  iteration  in  such  a  mixed  method  have  the  same  order  of  magnitude  as  the 
errors  in  the  second  iteration  of  the  extended  geometrical  method  proper. 

4)  If,  in  the  extended  geometrical  method  for  c»0.1  and  c»l,  the  second 
iteration  uses  only  second-order  partial  derivatives  (in  keeping  with  the 
pattern  of  the  first  Iteration),  the  errors  will  worsen  by  approximately  two 
orders  of  magnitude.  However,  this  does  not  lead  to  an  increase  in  the  number 
of  iterations  (namely  two)  resulting  from  the  current  cut-off  criterion. 


5)  In  this  example,  where  fourth-  and  higher-order  partial  derivatives  are 
zero,  a  procedure  utilizing  updated  second-order  partial  derivatives  in  the 
second  iteration  of  the  geometrical  method  agrees  perfectly  with  the  extended 
geometrical  method.  This  useful  verification  can  be  described  as  follows.  In 
general,  the  updated  second-order  partial  derivatives  can  be  expressed  by 


(du^) 

Par 
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since  here  the  quantities  represented  by  the  dots  are  zero,  the  equality  holds 
true  with  only  two  terms  on  the  right-hand  side.  But  this  is  precisely  the  form 
featured  by  the  extended  geometrical  method. 

6)  In  the  geometrical  and  the  extended  geometrical  methods,  as  well  as  in 
the  mixed  methods  of  the  comments  numbered  3  and  4,  the  solution  of  the  modified 
normal  equations  has  been  carried  out  using  both  the  Gauss  elimination  algorithm 
and  the  Choleskl  algorithm  for  positive-definite  matrices.  All  the  results. 
Including  the  Inverted  matrices  of  modified  normal  equations  in  all  iterations, 
have  been  found  identical  to  all  11  digits  printed.  This  illustrates  that  in 
the  usual  case  of  continuous  second-order  partial  derivatives  and  of  the 
solution  converging  to  a  minimum,  the  matrix  of  modified  normal  equations 
evaluated  in  the  neighborhood  of  the  final  point  P  is  positive-definite. 

The  above  example  illustrates  that  the  standard  method  converges  in  general 
only  for  small  residuals  (here  for  c=0.001  and  c=0.01).  When  larger  residuals 
are  present,  this  method  either  converges  very  slowly  or  diverges  (the  latter 
has  been  observed  for  c^O.l  and  c=l).  The  geometrical  method  is  particularly 
attractive.  It  is  relatively  simple,  yet  it  converges  in  two  iterations  for 
small  residuals  (c=0.001  and  c=0.01),  and  in  three  iterations  for  large 
residuals  (c=0.1  and  c=l).  In  the  case  of  large  residuals,  the  use  of  the 
extended  geometrical  method  further  reduces  the  number  of  iterations  to  two. 

This  reduction  has  been  accomplished  here  even  without  third-order  partial 
derivatives  in  the  second  iteration,  which  suggests  that  a  fourth  method  could 
be  adopted  (the  mixed  method  of  the  comment  4),  where  only  second-order  partial 
derivatives  evaluated  at  P  would  be  used,  whether  in  the  first  or  in  subsequent 
Iterations.  In  conclusion,  we  have  seen  that  the  standard  procedure  may  be  slow 
to  converge,  or  may  diverge  if  the  residuals  are  large,  whereas  the  use  of 
geometry  leads  to  fast  convergence,  represented  here  by  two  or  three  iterations 
depending  on  the  level  of  simplifications. 
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APPENDIX  B 


MODEL  SURFACE  AND  NOI»L  PLANE 

In  this  study,  the  underlying  space  encoapasslng  all  the  geoaetrlcal 
objects  Is  an  n-diaenslonal  flat  space  called  the  "observational  space".  The 


I* 

space  coordinate  systea  used  exclusively  Is  symbolized  by  {x  },  r>l ,2, . . . ,n,  In 
which  the  metric  tensor  Is 

e  9  constant  .  (B.l) 

*sr 


In  matrix  notation,  g  would  be  represented  by  a  positive-definite  matrix  of 

sr 

dimensions  (nxn).  Embedded  In  the  observational  space  Is  a  u-dlmenslonal 
surface  called  the  "model  surface".  The  surface  coordinate  system  Is  symbolized 

by  {u*},  a=l,2 . u.  The  surface  coordinates  are  considered  In  the  role  of  u 

Independent  variables,  without  any  specific  physical  meaning. 

In  the  observational  space,  the  model  surface  Is  defined  by  the  Gauss  form: 

x*"  =  x'*(u‘*)  ,  r  =  1,2 . n  ,  a  =  1,2, ...,u  .  (B.2) 

Equation  {B.2)  restricts  the  values  the  coordinates  x*^  can  have  In  order  to 
describe  points  on  the  model  surface.  This  form,  expressing  each  of  the  n  space 
coordinates  as  some  function  of  the  u  surface  coordinates.  Is  considered  known. 
It  Is  a  hlgher-dlmenslonal  analogue  of  the  Gauss  form  presented  in  Chapter  6  of 
[Hotlne,  1969],  describing  a  two-dimensional  surface  (u=2)  embedded  in  a  three- 
dimensional  space  (n>3). 

In  the  vicinity  of  an  "Initial"  point  P,  whose  model-surface  coordinates 
are  u®,  (B.2)  can  be  expanded  In  the  Taylor  series: 

x*"  -  x‘'(u®)  =  ^  <l/2)  (u“-u“)(u^-uj) 

+  (1/6)  (u®-u“)  (u^-u^)  (u’^-u^)  +  ...  ,  (B.3) 

where  0*'.,  ,  ...  are  coefficient  sets  evaluated  at  P.  The  set 

CL  Otp  m 

0  -  Is  symmetric  In  a  and  the  set  ®  -  Is  symmetric  In  a,  |3,  and  j.  etc. 
op  *pT 

The  differentiation  of  (B.3)  yields 

ax‘'/3u*  -  +  (1/2)  [d®{u^-uj)  +  (u®-u®)d^]  +  ...  . 

a  e  ap  6  0  o  « 

which,  upon  taking  advantage  of  the  symmetry,  becomes 
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3x'^/3u* 


•  *«  *  “is 

Similarly,  we  form 

a^x‘*/au“9u^  =  +  “'a^r +  ...  .  (b.s) 

a^x’'/au“3u^au’'  =  0**-  +  ...  .  (b.s) 

otpT 

etc.,  where  the  dots  Indicate  higher-order  terns. 

At  the  point  P,  where  u‘*=u^,  we  have  from  (B.3)-(B.6),  respectively: 


r  r 

X  =  X  : 

(B.3*) 

0 

-  r,-  o  .r 

ax  /au  =  A 

(X 

(B.4') 

a^x^/au^au^  =  , 

ap 

(B.5-) 

-3r,,a-(9-7'  ^r 

3  X  /au  au  au  =  ®  « 

otpT 

(B.8'  ) 

etc.,  confirming  the  coefficient  sets  in  the  Taylor  series  expansion  (B.3). 

Section  3.1  has  Introduced  as  the  "design  tensor",  expressible  via  the 
components  of  u  orthonormal  vectors  8,  J,  ...  spanning  the  hyperplane  tangent  to 
the  model  surface  at  P: 


-r,,a  .r  „r.  .r. 

ax/3u  =  A  =11  t  +  J  J  +  ..., 
a  a  <x 


(B.7) 


which  transforms  like  a  space  tensor  in  the  contravariant  indices,  and  like  a 
surface  tensor  in  the  covariant  indices.  Any  set  of  orthonormal  vectors 

r  oc 

spanning  this  hyperplane  is  acceptable.  In  [Hotlne,  1969]  the  tensor  ax  /du 


was  presented  for  n=3  and  u=2. 
tensors 


r  r 

On  the  other  hand,  n  ® 


are  not 


0(3’  aiSr’ 

For  example,  from  equation  6.14  in  [Hotlne,  1969]  it  follows  that 


a“,V3a«au“  -  a'  a»  .  a^  . 

r  r  S 

where  A  the  surface  covariant  derivative  of  A  with  respect  to  u  , 

o^graSr  “ 

tensor,  but  a  x  /au  au^-0  .  is  not,  due  to  the  Chrlstoffel  symbols. 

0/3 


is  a 


At  the  initial  point  P,  a  set  of  n  orthonormal  vectors  H ,  j .  u,  ... 

spans  the  n-dimenslonal  observational  space.  A  subset  of  u  vectors,  e,  J . 

spans  the  above  hyperplane,  while  the  remaining  n-u  vectors,  v .  are  normal 

to  it.  The  observational-space  metric  tensor  at  P  is  expressible  by 


4  1 


(B.8) 


sr 


e  +  j 
r  s 


Jr  " 


+  y  1/  +  . . .  . 

8  r 


The  model-surface  metric  tensor  at  P  is 


/3a  «  a  •’3  a 


A®  2  A** 

*sr 


which  follows  from  (B.7)  and  (B.8).  In  matrix  notation 
represented  by  a  positive-definite  matrix  of  dimensions 


a.  would  be 
pa 

(uxu) . 


{B.9) 


Next,  consider  another  u-dlmensional  surface  embedded  In  the  observational 

space  and  containing  the  point  P.  The  Independent  variables  {u*},  a=l,2 . u 

are  again  adopted  to  be  the  surface  coordinates,  and  their  values  at  P  are  again 
u^.  The  spatial  description  of  the  new  surface  is  provided  by 

x'*'  =  x'^(u®)  =  +  A**  (u^-u*)  .  (B.IO) 

o  a  o 

Contrasted  to  (B.4)-(B.6),  the  differentiation  now  yields 

3x'^/3u®  =  s  constant  ,  (B.ll) 

a 

3^x*'/3u^3u^  3  0  ,  3^x*'/3u*3u^3u^  =  0  ,  (B.12a,b) 

etc.  Similar  to  (B.3')  and  (B.4'),  at  P  we  have  x’*'»x*'  and  3x ' '^/3u*=A^, 

0  a 

r  r  r 

while  all  of  0  ,  ...  are  zero.  Since  A  in  (B.ll)  is  constant,  and 

ap  ap7  a 

g  Is  constant  by  definition,  the  metric  tensor  for  the  new  surface  Is 
3r 

=  Aa  g  A*'  s  constant  ,  (B.13) 

|3a  i3  ^sr  a 

Implying  that  the  surface  is  flat.  It  Is  thus  a  hyperplane,  which  we  call  the 
"model  plane".  Necessarily,  the  hyperplane  tangent  to  the  model  plane  at  P  is 
the  model  plane  itself. 

In  analogy  to  (B.7),  the  constant  tensor  In  (B.ll)  Is  written  in  the  form 

3x'^/3u®  =>  a''  =  6  ' +  j’'*  J'  +  ...  ,  (B.14) 

a  a  a 

where  the  orthonormal  vectors  t’,  J',  ...  lie  entirely  In  the  model  plane. 

Since  the  tensors  In  (B.7)  and  (B.14)  are  equal,  we  anticipate  that  the  vector 
sets  (,  J,  ...  and  t',  J',  ...  can  be  made  Identical,  and  thus  that  the  model 
plane  Is  tangent  to  the  model  surface  at  P.  This  follows  from  the  facts  below: 

1)  The  independent  variables  {u®}  represent  the  coordinate  system  for  both  the 
model  surface  and  the  model  plane; 
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2)  The  raodeJ -surface  metric  tensor  a.  at  P  is  the  same  as  Its  (constant) 

pa 

counterpart  for  the  model  plane;  and 

I* 

3)  The  model-surface  design  tensor  at  P  Is  the  same  as  Its  (constant) 
counterpart  for  the  model  plane. 


In  carrying  out  the  proof,  we  let  the  coordinate  differences  du*  represent 

an  Infinitesimal  model-plane  displacement  from  P  along  the  unit  vector  8*.  The 

contravariant  components  of  t  '  in  model-plane  coordinates  are  i'®*du**/ds,  where 

2  iS  (X 

ds  is  the  length  of  the  displacement  obtained  from  ds*=du^a-  du“.  Similarly, 
a 

one  formulates  J*  ,  etc.,  where 


«  '  fi  '  =  1 


° . 


Due  to  item  1,  the  set  du  can  also  represent  an  infinitesimal  displacement  in 
the  model  surface.  We  denote  the  length  of  this  displacement  by  di ,  and  a  unit 
vector  in  the  direction  of  the  displacement  by  ! .  The  contravariant  components 
of  1  in  model-surface  coordinates  are  2°=du*/d6.  Due  to  Item  2,  we  have  d2=ds 
and  thus  £^  =  C'®.  In  the  same  way.  .  Due  again  to  Item  2,  we  have 


£  =£',  J  »J’ . so  that 

a  a  a  a 


£  =  1  . 
a 


J. 


0  ,  . 


Accordingly,  the  thus  constructed  vectors  «,  J,  ...  behave  like  model-surface 
orthonormal  vectors  at  P. 


With  regard  to  the  spatial  configuration  of  £.  j .  Item  3  reads 

a  •'a  a 

where  (B.7)  and  (B.i4)  have  been  utilized  together  with  the  fact  that  = 

j'  =  j .  Upon  successive  contractions  with  £**,  J* . It  follows  that 

a  a 

r  r  r  r 

£  -  £  ’  ,  J  -  J’  .  (B.15) 

Prom  the  equality  of  their  space  components  we  deduce  that  the  vectors  £  and  £', 
J  and  j',  etc.,  are  Identical.  Since  £,  J.  ...  are  tangent  to  the  model  surface 
at  P,  and  £’,  j',  ...  He  entirely  In  the  model  plane,  the  model  plane  defined 
by  (B.IO)  is  seen  to  be  tangent  to  the  model  surface  defined  by  (B.3).  The 
point  of  tangency,  P,  is  described  by  u®  In  the  common  surface  coordinate 
system  (u®),  and  by  x^  =  x*'(u®)  in  the  space  coordinate  system  (x'’). 
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APPENDIX  C 


LINKS  BETNEEN  OBSERVATI(MIAL-SPACE .  NOOEL-PLANB. 
AND  NODEL-SORFACE  CO(XiDINATE  DIFFERENCES 


If,  as  in  the  current  study,  the  underlying  observational  space  is  flat  and 

r* 

the  space  coordinate  systea  {x  }  is  such  that 


g  3  constant  , 
sr 

in  analogy  to  Cartesian  coordinates  it  can  be  shown  that 

Ax'*  =  ,  (C.l) 

r 

where  the  set  Ax  contains  the  observational-space  coordinate  differences 
between  two  points,  and  Ap  symbolizes  the  difference,  in  contravarlant 
components,  between  the  position  vectors  of  these  two  points.  The  set  Ap^  gives 
rise  to  the  tensor  equation 


Ap**  =  As  k**  ,  (C.2) 

I* 

where  k  represents  the  contravarlant  coaponents  of  the  unit  vector  k  in  the 
direction  of  Ap,  and  As  is  the  length  of  Ap.  Just  as  in  Cartesian  coordinates, 
the  position  vector  p  belonging  to  any  point  in  space  can  be  freely  parallel- 
transported  from  that  point  to  an  arbitrary  location,  and  the  same  can  be  said 
about  Ap  and  k. 

Since  we  shall  use  exclusively  {x**}  as  the  observational-space  coordinate 
system,  we  write  Ax^  for  Ap^ .  although  in  general  coordinates  (C.l)  would  not 
hold  true  (it  is  not  a  tensor  equation).  Prom  (C.2)  we  deduce 


As  =  Ax'"  k  =  (l/As)  Ax*'  Ax 
r  r 

or 

»2  ,3.  ,3  .r 

A8>AxAx=-Axg  Ax, 
8  ‘sr 

where 


(C.3) 


Ax  >  g  Ax  . 
s  sr 

Equation  (C.4)  is  written  in  lieu  of  the  tensor  equation  Ap^-g^^Ap**. 


(C.4) 
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SiaiJar  considerations  apply  for  the  model  plane  with  a-  aconstant  as  given 

pa 

by  (B.13).  In  analogy  to  (C.l)  and  (C.2),  we  thus  write 


Au 


A  M 


As  k“ 


(C.5) 

(C.6) 


Here  the  symbols  Ap.  As.  and  k  are  used  Independently  of  the  same  symbols 
appearing  In  (C.1,2).  In  analogy  to  (C.4),  we  also  have 


''“e  ■  “s. 


representing  the  tensor  equation  Ap. =a-  Ap  . 

P  pa 


(C.7) 


In  linking  observational-space  and  model-surface  coordinate  differences  we 
shall  proceed  via  the  model  plane,  where  we  can  take  advantage  of  the  property 
(C.5).  Let  Ax'  reptesent  a  model-plane  vector  emanating  from  P.  which,  In 
Fig.  1,  is  the  vector  symbolized  by  the  arrow  PP^.  Upon  considering  i,  j,  ... 
to  be  orthonormal  vectors  spanning  the  model  plane,  it  follows  that 

Ax'''  =  ae''  +  bj''+...  .  (C.8) 


where  a,  b,  ...  are  scalars.  If  the  same  vector  should  be  described  by  model- 
plane  rather  than  observational-space  components,  we  denote  It  Au  and  write 

Au  »a8*-bj+....  (C.9) 


Due  to  (B.14),  presented  here  as 

a'”  =  e''  e  t  J  +  .  .  .  .  (C.IO) 

a  a  a 

one  obtains 

Ax'^'  =  Au“  .  (C.ll) 

where  Ax'  symbolizes  the  contravarlant  components  of  the  above  vector  as  well 
as  the  corresponding  coordinate  differences  in  the  observational-space 
coordinates,  and  Au*  symbolizes  such  components  and  such  coordinate  differences 
in  the  model-plane  coordinates. 

Since  both  (C.8)  and  (C.9)  hold  true  also  when  the  Indices  are  lowered,  we 
readily  deduce  that 

AUp  =  a|  g^^  Ax'**  .  (C.12) 
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This  relation  can  also  be  obtained  froa  (C.ll),  upon  contracting  the  latter  with 

A^g  and  recalling  (C.7)  together  with  (B.13).  As  a  special  benefit,  if  Ax' 
p  sr 

were  replaced  by  a  general  space  vector  Ax.  the  left-hand  side  of  (C.12)  would 
describe  a  projection  of  that  vector  onto  the  aodel  plane. 

To  conflra  this  property,  suppose  that  the  contravariant  coaponents  of  a 
r  r 

space  vector  are  Ax'  +Ax"  ,  where 
Ax"**  =•  q  w*"  +  .  . .  . 

In  terns  of  Fig.  1,  we  can  Imagine  Ax"  as  a  vector  perpendicular  to  the  model 
plane  at  P^.  It  then  follows  that 

%  =*  A®  ggj.  (Ax'^'  +  Ax""')  =  a|  gg^  Ax'^  .  (C.13) 

In  other  words.  If  the  components  Ax'^<-Ax"*^  are  known,  (C.13)  yields  the  model- 
plane  covariant  components  of  the  projected  vector,  Au*.  The  observational- 

”  p 

space  contravariant  components  of  the  projected  vector,  Ax'  ,  can  be  obtained 
via  (C.ll)  upon  first  utilizing 

Au*  -  Au .  ,  (C . 14 ) 

p 

where  a*^  Is  the  aodel-plane  associated  metric  tensor  (such  that  a**^ap^=fi®). 

One  feature  of  equations  (C.13)  and  (C.14)  is  that  they  demonstrate  how  the 
observational-space  components  of  a  vector  can  be  converted  in  ■>  model-plane 
coordinate  differences. 


We  shall  next  present  a  spatial  link  between  the  model  surface  and  the 
model  plane  by  relating  corresponding  points  in  the  two  surfaces.  In  Fig.  1, 
a  pair  of  such  points  Is  depicted  by  P  and  P^.  Their  correspondence  is 
understood  In  the  sense  that  the  same  set  of  values  u“  represents  the  model - 
surface  coordinates  of  P  and  the  model-plane  coordinates  of  P^.  In  the 
observational  space,  the  points  P  and  P^  are  identified  by  their  position- 
vector  components  x*"  and  x'^,  respectively.  By  subtracting  from  either  of  the 

J* 

latter  the  position-vector  components  of  P,  x^,  we  have 


r  r 

X  -  X 

o 


These  sets  of  contravariant  components,  which  represent  also  sets  of  coordinate 
differences,  identify  the  vectors  Ax  and  Ax'  symbolized  In  Pig.  1  by  the  arrows 
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PP  and  PP^.  respectively .  Upon  denoting 


the  coflipatison  of  (B.3)  and  (B.IO)  reveals  that 


r  r  r 

Ax  «  Ax*  +  V 

where  the  contravarlant  components  v  of  the  vector  v  are  given  as 
v*"  =  (1/2)  Au®  +  (1/8)  Au“  Au^  Au’'  +  .  .  .  . 

The  vector  v  provides  the  desired  link  between  the  two  surfaces. 


(C.15) 

(C.16) 


In  considering  that  the  symbols  x  are  Interchangeable  here  with  p,  the 

r  r 

relation  (C.15)  Is  a  tensor  equation.  On  the  other  hand,  since  0  ®  , 

«p  <xpr 

...  are  not  tensors,  (C.16)  Is  not  a  tensor  equation,  but  this  does  not  detract 
from  Its  usefulness.  We  observe  that  tensor  considerations  have  been  Important 
mainly  in  deriving  (C.13),  where  the  symbols  Au  and  Ax'  represent  the  same 
geometrical  object.  The  tensor  contraction  with  a®^  In  (C.14)  then  yields  Au®. 
which  Identifies  both  the  set  of  contravarlant  components  of  the  vector  PP^  and 
the  set  of  model-plane  coordinate  differences  between  Pg  and  P.  Due  to  the 
definition  of  (u®),  the  set  Au®  represents  also  model-surface  coordinate 
differences  between  P  and  P. 


APPENDIX  D 


MATRIX  NOTATION 

If  needed,  the  results  of  this  study  can  be  transcribed  into  Matrix 
notation.  This  Is  greatly  facilitated  by  the  close  correspondence  between 

p 

tensor  contractions  and  matrix  multiplications.  In  the  cases  of  0  .  and 
r 

«  a  .  respectively,  we  introduce  three-  and  four-dinensional  arrays  into  the 
ctpT 

matrix  context.  In  transcribing  the  corresponding  quantities,  the  indices  will 

be  dropped.  The  letter  symbols  will  be  retained,  but  *  will  be  attached  to 

1*8 

quantities  corresponding  to  purely  covariant  tensors.  Thus,  the  tensors  g 

a®^,  dx"^,  Ax^.  Au®,  etc.,  will  become  g.  A,  a,  dx.  Ax,  Au,  etc.,  while  the 
a 

tensors  g  .  a.  ,  AUa,  etc.,  will  become  g* ,  a*,  Au*,  etc.  The  first  index 
Identifies  the  rows  and  the  second  Index,  if  present,  identifies  the  columns. 

If  the  Indices  are  mixed,  as  in  A^,  the  upper  (contravariant)  index  io 
regarded  as  the  first  and  the  lower  (covariant)  index  is  regarded  as  the  second, 
r  r 

The  sets  0  .  and  <D  .  will  be  transcribed  respectively  as  a  three- 
ap  apT 

dimensional  array  0  of  dimensions  (nxuxu)  and  a  four-dimensional  array  ®  of 
dimensions  (nxuxuxu).  In  analogy  to  the  matrix  A  of  dimensions  (nxu)  evaluated 
at  P  according  to 

A  s  [ax/au^  ax/au^,  . . .  ]  . 

1  2 

where  x  represents  the  observables  and  u  ,  u  ....  are  the  individual  parameters, 
the  arrays  0  and  0  are  formed  at  P  as 

n  2  [aA/au^  aA/au^,  ...  ]  .  0  =  [an/au^  an/au^,  ...  ]  . 

Upon  considering  the  order  of  contractions,  it  becomes  clear  that  the  expression 
9  T 

A^g^r  corresponds  to  A  g*  (and  not.  for  example,  to  Ag*  which  is  not  even 

defined).  The  contractions  of  the  kind  D^-(Au^)  and  {5x  )®a  will  be 

j  0(3 '  '  '  3 '  0or 

transcribed  as  fl(Au)  and  (dx*)  w,  respectively. 

r  r 

In  dealing  with  quantities  related  to  the  observables,  such  as  x^.  x  . 
r  r  r  “T 

x^,  dx  ,  Ax  .  dx  ,  etc.,  which  in  the  geometrical  context  of  the 
observational  space  represented  sets  of  coordinates  or  sets  of  coordinate 
differences  equivalent  to  tensors,  we  again  drop  the  Indices  and  use  the  symbols 
Xq,  X.  x^,  dx.  Ax.  dx,  etc.  The  column  vectors  (of  n  elements)  x^,  x,  and  x^ 
denote  respectively  the  actual  observations,  the  adjusted  observations,  and  the 
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observables  consistent  with  the  Initial  set  of  parameters.  The  column  vectors 
Sx.  Ax,  and  <5x,  the  latter  containing  the  negative  residuals,  follow 
respectively  from  the  relation  below  (8),  from  (4),  and  from  (14): 

dx  2  x„  -  X  : 

Q  o 

Ax  =  X  -  X  .  dx  3  x„  -  X  =  <5x  -  Ax  . 

o  Q 

Quantities  related  to  the  parameters,  such  as  u®,  u®,  Au®,  etc.,  which  in  the 
geometrical  context  of  the  model  plane  are  sets  of  coordinates  or  set*:  of 
coordinate  differences  equivalent  to  tensors,  are  symbolized  by  u,  u^.  Au.  etc. 
The  column  vectors  (of  u  elements)  u  and  denote  respectively  the  adjusted 
parameters  and  the  initial  values  of  parameters,  while  the  column  vector  Au 
follows  from  (lb)  as 

Au  =•  u  -  u 

o 

Introducing  also  quantities  in  parentheses,  in  analogy  to  Ax  and  Sx  we 

have 

(Ax)  =  (x)  -  x^  ,  (<5x)  =  Xq  -  (x)  =  5x  -  (Ax)  . 

Moreover,  in  analogy  to  Au  we  write 
(Au)  =  (u)  -  u^  . 

The  intermediate  values  such  as  (u)  and  (x)  follow  from  these  equations  as 
(u)=u  +{Au)  and  (x)=x  +(Ax).  Similar  transcriptions  apply  also  for  the 
relationships  involving  A(Au  )  and  A(Ax  ),  such  as  presented  in  (15)  and  (19a). 
respectively,  but  it  is  unnecessary  to  write  them  explicitly  since  the  pattern 
explained  above  is  general  and  straightforward. 

We  next  transcribe  the  key  formulas  associated  with  the  initial  iteration 
of  Section  3.1.  Equations  (6),  (8),  and  (9)  correspond  to 

a*  =  g*  A  ,  (D.l) 

(Au*)  =  A^  g*  <5x  ,  (Au)  »  a  (Au*)  .  (D.2a,b) 

where  a  is  the  Inverse  of  a*.  Further,  the  formulas  (10a, b)  become 

g'  »  A  a  A^  ,  g*'  =  g»  g'  g*  .  (D.3a.b) 
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while  the  relations  corresponding  to  (11a, b)  are 

g"  =  g  -  g'  .  g*"  =«  g*  -  g*'  .  (D.4a.b) 

The  adjustment  role  of  the  matrices  g* ,  g*',  g",  and  g*"  of  dimensions  (nxn)  has 
already  been  explained.  Similar  to  the  relationship  between  the  matrices  a  and 
a*  of  dimensions  (uxu),  the  matrices  g  and  g*  of  dimensions  (nxn)  are  Inverses 
of  each  other. 

Finally,  we  present  the  system  (23a~c)  In  matrix  notation: 

[(a*)  -  (<5x*)^  0  -  (<5x*)”^  ®  (Au)]  A(Au)  =  (A)’*  (Ax*)  ,  (D.5a) 

where 

(a*)  =  (A)^  g*  (A)  .  ((5x*)  =  g*  (dx)  .  (D.Sb.c) 

We  call  (D.5a)  the  modified  normal  equations,  and  the  matrix  on  the  left-hand 
side,  denoted  (a*),  the  matrix  of  modified  normal  equations.  For  a  repeated 
application  of  (D.5a-c),  (Au)  is  updated  through  the  algebraic  addition  of 
A(Au).  All  the  other  quantities  in  parentheses  are  updated  as  well,  while  the 
arrays  (1  and  <6  are  treated  as  constant.  Depending  on  the  assumption  about  the 
size  of  the  residuals.  (D.5a)  can  be  simplified  by  leaving  out  the  term  with  ®. 
If  both  terms  wltii  fi  or  ®  are  left  out,  the  resulting  equation  characterizes  the 
standard  adjustment  approach.  The  final  updated  values  are  denoted  by  overbars. 
The  variance-covariance  and  the  weight  matrices  of  adjusted  quantities  are 
represented  by  (D.l),  {D.3a,b).  and  (D.4a.b),  with  all  the  symbols  except  g  and 
g*  overbarred,  and  with  a  computed  as  the  inverse  of  a*. 
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